-------- Original Message -------- Subject: repost of an older post concernic custom R functions Date: Wed, 8 Sep 2010 13:33:42 +0100 From: Milos Blagojevic <[email protected]> To: <[email protected]> Hello to all morphmet users, This is just a repost of my email from a couple of months back which I never saw reach mailing list. I don`t know if this is due to my email account or maybe to the list itself. Either way I am sorry if You already read this one. Recently, I have been working on some custom functions in R (on, I guess introductory level) in order to perform bootstrap analysis of Procrustes distances. Firstly i wrote this function to calculate partial, Procrustes and full Procrustes distances between two coordinate matrices (and a rewriting of it to accept vectors instead of matrices): procDist <- function(M1,M2) { num <- dim(M1)[1] diffe1 <- rep(NA,num) diffe2 <- rep(NA,num) for(i in 1:num) {diffe1[i] <- (M1[i,1]-M2[i,1])^2 diffe2[i] <- (M1[i,2]-M2[i,2])^2} diffe <- c(diffe1,diffe2) Dp <- sqrt(sum(diffe)) Pdrho <- 2*asin(Dp/2) Df <- sin(Pdrho) list(Dp = Dp, Pd = Pdrho, Df = Df) } procDistvect <- function(a,b) { num <- length(a) diffe <- rep(NA,num) for(i in 1:num) {diffe[i] <- (a[i]-b[i])^2} Dp <- sqrt(sum(diffe)) Pdrho <- 2*asin(Dp/2) Df <- sin(Pdrho) list(Dp = Dp, Pd = Pdrho, Df = Df) } These both work fine...but now I wanted to perform a bootstrap test to estimate significance of Procrustes distance between mean male and female shape in my sample. y1 <- c(1:50) mshapesMA <- array(NA,dim=c(19,2,4000)) - empty arrays to write calculated bootstrapped mean shapes mshapesFE <- array(NA,dim=c(19,2,4000)) for(i in 1:4000) { resamples <- lapply(1:4000, function(i) sample(y1, replace=T)) -these are the resamples (I have 50 individuals in a sample) mshapesMA[,,i] <- mshape(procCCall[,,resamples[[i]]][,,1:39]) - calculations for mean shapes of bootstraps (39 males and 11 females) mshapesFE[,,i] <- mshape(procCCall[,,resamples[[i]]][,,40:50])} bootPD <- rep(NA,4000) for(i in 1:4000) {bootPD[i] <- procDist(mshapesMA[,,i],mshapesFE[,,i])$Df} - taking calculated Procrustes distances the results I get in bootPD variable seem fine at first and are normally distributed. But there is a flaw...all calculated PDs tend to be smaller than the real, one-by-one calculated PDs. I have no slightest idea why this happens but i would be grateful for any help in solving this problem. Also it is possible to use the famous boot R library, but the problem here is my inability to adjust procDist function to work with boot function. This boot accepts only such functions that are built with two inputs, one the matrix of coordinates and the other the row index that boot is shuffling. This version: procDistvectboot <- function(a,d) { num <- length(a[1,]) diffe <- rep(NA,num) for(i in 1:num) {diffe[i] <- (a[d,i]-a[d,i])^2} Dp <- sqrt(sum(diffe)) Pdrho <- 2*asin(Dp/2) Df <- sin(Pdrho) list(Dp = Dp, Pd = Pdrho, Df = Df) } just doesn`t work. I can`t find a solution to indexing(d) of coordinate (XXXXYYYY, rows are numbers of individuals, and columns are landmarks) matrix commonly used in R. I always find a conflict between row and column number indexing in loops as well as in all apply-type code. The key is to calculate the differences between corresponding landmarks in first two rows of (a) matrix, and then between the second and third row and so on. Indexing of (a) always fails and doesn`t seem to work in any of the configurations that I could remember. This bothers me so much lately and help is greatly needed, Best regards Milos Blagojevic,PhD student Faculty of science Depatment for biology and ecology Radoja Domanovica 12 34000 Kragujevac Serbia e-mail:[email protected] - morphmet [email protected]
