[R] Simplify formula for iterative programming

[Christoph Buser]

Dear Stef

Please check my code below carefully, because it is Friday
afternoon and there might be mistakes in it:

## initialize two random vectors (length n = 100)
a <- rnorm(100)
b <- rnorm(100)
## computation of your formula (see ?outer which is very
## useful here) 
(c <- sum(sqrt(outer(a,a,"-")^2 + outer(b,b,"-")^2)))
## expand a and b to length n+1
a1 <- c(a,rnorm(1))
b1 <- c(b,rnorm(1))
## computation calculating the whole sum
(c1 <- sum(sqrt(outer(a1,a1,"-")^2 + outer(b1,b1,"-")^2)))
## computation using the result of the smaller vectors
c + 2*sum(sqrt((a1 - a1[length(a1)])^2 + (b1 - b1[length(b1)])^2))

Regards,

Christoph Buser

--------------------------------------------------------------
Christoph Buser <buser at stat.math.ethz.ch>
Seminar fuer Statistik, LEO C13
ETH (Federal Inst. Technology)	8092 Zurich	 SWITZERLAND
phone: x-41-44-632-4673		fax: 632-1228
http://stat.ethz.ch/~buser/
--------------------------------------------------------------


Stefaan Lhermitte writes:
 > Dear R-ians,
 > 
 > I am looking for the simplification of a formula to improve the 
 > calculation speed of my program. Therefore I want to simplify the 
 > following formula:
 > 
 > H = Si (Sj ( sqrt [ (Ai - Aj)² + (Bi - Bj)²  ] ) )
 > 
 > where:
 > A, B = two vectors (with numerical data) of length n
 > sqrt = square root
 > Si = summation over i  (= 0 to n)
 > Sj = summation over j (= 0 to n)
 > Ai = element of A with index i
 > Aj = element of A with index j
 > Bi = element of B with index i
 > Bj = element of B with index j
 > 
 > n is not fixed, but it changes with every run for my program. Therefore 
 > for I am looking for a simplication of h in order to calculate it when 
 > my A and B get extendend by 1 element (n = n + 1).
 > 
 > I know a computional simplified formula exists for the standard 
 > deviation (sd) that is much easier in iterative programming. Therefore I 
 > wondered I anybody knew about analog simplifications to simplify H:
 > 
 > sd = sqrt [ ( Si (Xi - mean(X) )² ) /n  ]  -> simplified computation -> 
 > sqrt [ (n * Si( X² ) - ( Si( X ) )² )/ n² ]
 > 
 > This simplied formula is much easier in iterative programming, since I 
 > don't have to keep every element of X.
 > E.g.: I have a vector X[1:10]  and I already have caculated Si( X[1:10]² 
 > ) (I will call this A) and Si( X ) (I will call this B).
 > When X gets extendend by 1 element (eg. X[11]) it easy fairly simple to 
 > calculate sd(X[1:11]) without having to reuse the elements of X[1:10].
 > I just have to calculate:
 > 
 > sd = sqrt [ (n * (A + X[11]²) - (A + X[11]²)² ) / n² ]
 > 
 > This is fairly easy in an iterative process, since before we continue 
 > with the next step we set:
 > A = (A + X[11]²)
 > B = (B + X[11])
 > 
 > Can anybody help me to do something comparable for H? Any other help to 
 > calculate H easily in an iterative process is also welcome!
 > 
 > Kind regards,
 > Stef
 > 
 > ______________________________________________
 > R-help at stat.math.ethz.ch mailing list
 > https://stat.ethz.ch/mailman/listinfo/r-help
 > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html