[R] help in matlab - r code

Susanne Schmidt suse_semperula at yahoo.de
Mon Mar 29 17:07:08 CEST 2010

Dear Marta,

I did it in Matlab, and fiddled around with R code until I had *almost* the same result. The "almost" is probably due to R handling the picture values (ranging from 0 to 1) differently than Matlab (ranging from 0 to 255), and simply multiplying the R picture values by 255 did NOT result in exactly the same values as the Matlab values. [what seems white in the picture is 245 in Matlab, although values potentially range to 255, and white is 0.9642549 in R, which multiplied by 255 gives 245.12, e.g.]
But maybe the precision of this solution is good enough for you ..

The corr2 demand from matlab is a 2D correlation coefficient - the R command cor works elementwise, and is not the solution here.
Below I tried to implement the formula given in the following matlab page:

Maybe somebody on the list has a nice idea how to make the code more elegant

This is the complete code in R

setwd("D:/   wherever ")
x <- read.jpeg("  whichever  .jpg") #open image

plot(x) #plot image
x <- rgb2grey(x) #convert to greyscale
plot(x) # check ;-) the image is in grey scale

ImageWidth = dim(x)[2] #number of col
MaxOffset = 99; #defined variable
ImageWidthToProcess = ImageWidth-MaxOffset; #col-defined variable

## this one does NOT work because matrices not square:
for(k in 1: MaxOffset) {
   OffsetPlaquette <- x[  , c((1+ k)   :  (ImageWidthToProcess + k))]
   dataToProcess <- x[,c(1:ImageWidthToProcess)]
 AutoCData[k] <-  mantel(OffsetPlaquette, dataToProcess)
## END this one does not work because matrices not square

AutoCData <- rep(0, MaxOffset)
sumBothM <- rep(0, MaxOffset)
sum1stMsq <- rep(0, MaxOffset)
sum2ndMsq <- rep(0, MaxOffset)
for(k in 1: MaxOffset) {
   OffsetPlaquette <- x[  ,(1+k) : (ImageWidthToProcess + k)]
   dataToProcess <- x[,c(1:ImageWidthToProcess)]
   meanM <- mean(OffsetPlaquette); meanM2 <- mean(dataToProcess)
   for(j in 1:dim(dataToProcess)[2]){
     for(i in 1:dim(OffsetPlaquette)[1]){
        sumBothM[k]  <- sumBothM[k]  +  (OffsetPlaquette[i,j]-meanM)*(dataToProcess[i,j]-meanM2)
        sum1stMsq[k] <- sum1stMsq[k] +  (OffsetPlaquette[i,j]-meanM)2
        sum2ndMsq[k] <- sum2ndMsq[k] +  (dataToProcess[i,j]-meanM2)2
 AutoCData[k] <-  sumBothM[k]/(sqrt(sum1stMsq[k] *  sum2ndMsq[k]))


Best wishes,

Do You Yahoo!?
Sie sind Spam leid? Yahoo! Mail verfügt über einen herausragenden Schutz gegen Massenmails. 

More information about the R-help mailing list