Good day! I need to implement one of Knut's tests. Chose the criterion of serial correlation. The second volume of Knut gives the formula

formula

I'm trying to calculate this correlation coefficient from it and I get it at N = 50 equal to 1.04. But it cannot be, by module, more than one! ... Help, please, sort it out ...

A piece of code:

const int M=N; unsigned int null_num, first_num, numerator_n, numerator_2, denominator_n, denominator_2; long double C; null_num = (a*x+c)%4294967296; numerator_n=0; numerator_2 = null_num; denominator_n = null_num*null_num; for (i=0; i<M; i++) { first_num = (a*null_num+c)%4294967296; numerator_n = numerator_n + null_num*first_num; numerator_2 = numerator_2 + first_num; denominator_n = denominator_n + first_num*first_num; null_num = first_num; } numerator_n = numerator_n + null_num*((a*x+c)%4294967296); numerator_2 = numerator_2 + null_num; denominator_2 = numerator_2; denominator_n = denominator_n + null_num*null_num; C = (double)((numerator_n*M)-(numerator_2*numerator_2))/(double)((denominator_n*M)-(denominator_2*denominator_2)); printf("Coefficient: %f \n\n", C);
  • what language? from 4294967296? - Gedweb
  • 4294967296 - 2 ^ 32 four bytes, 2 ^ (8 * sizeof (int)) - zb '11
  • and where is he from in the formula - VladD
  • it's only the author knows, there are no divisions in the formula with the remainder, I think some kind of overflow protection but I am not sichnik - zb '11
  • @ 777Julia777 values ​​and types a, x, c provide - zb '11

1 answer 1

Your test is incorrectly implemented.

Do something like this:

 double sumUV = 0, sumU = 0, sumV = 0, sumU2 = 0, sum V2 = 0; for (int i = 0; i < n; i++) { double u = GetNextU(i); double v = GetNextV(i); sumU += u; sumV += v; sumU2 += u*u; sumV2 += v*v; sumUV += u*v; } double coeff = (n * sumUV - sumU * sumV) / sqrt((n * sumU2 - sumU * sumU) * (n * sumV2 - sumV * sumV));

This, unfortunately, is a rather naive code: due to inaccuracies in rounding, the result can be catastrophically inaccurate. Read more in Knut, 4.2.2. The example after formula (15) shows that when calculating the root in the denominator, you may encounter a negative number, despite the fact that it is mathematically impossible.

Update With the new formula you have V_i = U_{i+1} , you consider the sequence correlation to yourself. The code can be simplified.

 double sumU = 0, sumU2 = 0, sumUShifted = 0; double firstU = GetNextU(0); double prevU = firstU; for (int i = 0; i < n; i++) { double nextU = i == (n - 1) ? firstU : GetNextU(i + 1); sumU += prevU; sumU2 += prevU*prevU; sumUShifted += prevU*nextU; prevU = nextU; } double coeff = (n * sumUShifted - sumU * sumU) / (n * sumU2 - sumU * sumU);
  • Well, from and not what 4294967296 - Gedweb
  • @GedWeb: 4294967296 goes into the implementation of GetNextU , but there, even speaking, it is possible without it. - VladD
  • one
    @VladD, and what's the point in auto? All the same, because everything is further in double is considered. You can require GetNextU () to return double. But, IMHO it will not save the vehicle from loss / overflow at N = 50. No wonder a decent part of the second volume is given to calculations in rational numbers. - avp
  • @avp: you're right, I'll clean it up now. In TS, they were integers, but if you start a variable as a double, it will be more correct. - VladD
  • 2
    @ 777Julia777: great, good luck to you! But beware of accumulating inaccuracies. - VladD