[R-sig-eco] Question about a resampling procedure with LME

[Thiago Gonçalves-Souza]

Hello all,



First, I'm sorry for this long question, guys.


I’d like to know if the diversity of waterbirds affect zooplankton
diversity in natural temporary pools over two seasonal periods (wet and
dry). I’ve sampled zooplankton in 38 10x15 cm plots on the pool’s shore
sediment (in the same place where birds usually) in both seasonal periods.
I used a model of bird feet and simulated them walking on the sediment of
each plot. Then, I washed the bird feet (of each plot separately) into
containers with filtered water. The zooplankton that emerged from this
water were identified.


Unfortunately, I was not able to collect birds at the same plot (and at the
same time) that I simulated their footprint. Rather, I’ve made 38 samples
(during 9 months immediately after each seasonal period) of bird
communities (direct observation on each pool).


Thus, I do not believe that can use any classic model (such as anova
designs) because birds and zooplankton were not collected at the same
sampling unit.

To do this, I thought to randomly select *10* plots of each combination of
zooplankton sampling, that is: 10 samples from pool 1/wet season, 10 from
pool 1/dry season, 10 from pool 2/wet season and 10 from pool 2/dry season.
I did the same with bird samples, selecting 10 samples by pool and season
combination. Then, we run a lme model using species richness of birds
(independent variable) and zooplankton (dependent variable), and season as
independent variable. The variables were nested in pool. I repeated (for
loop) this procedure 1000 times and calculated the number of times that the
model was significant (p < 0.05) and the accompanying F values (for bird
richness and season effects).


I’d like to hear from you if this procedure make any sense. Please, see
below a reproducible example. Thank you very much!


Note: I’m pretty sure that the for loop could be shorter ;)



library(nlme)



#### Fake da ta



## zooplankton table

ric.zoo <- rnorm(120, 20, 2)

pool <- factor(rep(c("pool1", "pool2"), 60))

season <- factor(rep(rep(c("dry", "wet"), each=30), 2))

tab <- data.frame(pool, season, ric.zoo)



## birds table

ric.bird <- rnorm(76, 10, 2)

pool_b <- factor(rep(c("pool1", "pool2"), 38))

season_b <- factor(rep(rep(c("dry", "wet"), each=19), 2))

tab_b <- data.frame(pool=pool_b, season=season_b, ric.bird)



## Dataset to be used in for loops

bird.p.val <- numeric()

season.p.val <- numeric()

bird.f <- numeric()

season.f <- numeric()

pool.sort <- rep(c("pool1", "pool2"), each=20)

season.sort <- rep(rep(c("wet", "dry"), each=10), 2)

factor.sort <- data.frame(pool=pool.sort, season=season.sort)

nreps <- 1000





for(i in 1:nreps){

  subset1 <- subset(tab, pool=="pool1"&season=="wet")

  subsample1 <- subset1[sample(nrow(subset1), 10),][,3]

  subset2 <- subset(tab, pool=="pool1"&season=="dry")

  subsample2 <- subset2[sample(nrow(subset2), 10),][,3]

  subset3 <- subset(tab, pool=="pool2"&season=="wet")

  subsample3 <- subset3[sample(nrow(subset3), 10),][,3]

  subset4 <- subset(tab, pool=="pool2"&season=="dry")

  subsample4 <- subset4[sample(nrow(subset4), 10),][,3]

  subbird1 <- subset(tab_b, pool=="pool1"&season=="wet")

  subsbird1 <- subbird1[sample(nrow(subbird1), 10),][,3]

  subbird2 <- subset(tab_b, pool=="pool1"&season=="dry")

  subsbird2 <- subbird2[sample(nrow(subbird2), 10),][,3]

  subbird3 <- subset(tab_b, pool=="pool2"&season=="wet")

  subsbird3 <- subbird3[sample(nrow(subbird3), 10),][,3]

  subbird4 <- subset(tab_b, pool=="pool2"&season=="dry")

  subsbird4 <- subbird4[sample(nrow(subbird4), 10),][,3]

  ric.sort.zoo <- c(subsample1, subsample2, subsample3, subsample4)

  ric.sort.birds <- c(subsbird1, subsbird2, subsbird3, subsbird4)

  df40 <- data.frame(ric.zoo=ric.sort.zoo, ric.bird=ric.sort.birds,
factor.sort)

  lm.nested <- lme(ric.zoo~ric.bird+season, random=~1|pool, data=df40,
method="REML")#MODELO

  mod <- anova.lme(lm.nested, type="sequential", adjustSigma = FALSE)

  bird.p.val[i] <- mod$`p-value`[2]

  season.p.val[i] <- mod$`p-value`[3]

  bird.f[i] <- mod$`F-value`[2]

  season.f[i] <- mod$`F-value`[3]

}





sum(bird.p.val<0.05 )/nreps

sum(season.p.val<0.05)/nreps

-- 

*Thiago Gonçalves-SouzaProfessor Adjunto*

*Universidade Federal Rural de Pernambuco (UFRPE)*
*Departamento de Biologia / Área de Ecologia *

*E-mail alternativo: thiagoaracno at gmail.com <thiagoaracno at gmail.com>Home
page:<http://thiagocalvesouza.wix.com/ecoffun
<http://thiagocalvesouza.wix.com/ecoffun>**>*
*Google Scholar:
<http://scholar.google.com.br/citations?user=TjaP2l8AAAAJ&hl=pt-BR
<http://scholar.google.com.br/citations?user=TjaP2l8AAAAJ&hl=pt-BR>>*

	[[alternative HTML version deleted]]