[R-sig-eco] Comparing multiple distributions within treatments

Law, Jason Jason.Law at portlandoregon.gov
Thu Apr 23 01:29:15 CEST 2015


Jeff,

You can treat pair as a fixed effect and analyze the data using Poisson glm:

d <- gather(my.data2, 'tolerance', 'count', num_range('t', 1:4))
summary(g <- glm(count ~ (my.pair * my.treat) + tolerance:(my.pair + my.treat), data = d))

See section 7.3 of MASS.

The same section discusses a proportional odds approach to a similar problem. Something like:

summary(polr(tolerance ~ my.pair + my.treat, data = d, weights = count))

You'd want to check the proportional odds assumption.

There seems to be an 'ordinal' package which seems like it's right up this ally. 

Jason

-----Original Message-----
From: R-sig-ecology [mailto:r-sig-ecology-bounces at r-project.org] On Behalf Of Holland, Jeffrey D
Sent: Wednesday, April 22, 2015 1:30 PM
To: r-sig-ecology at r-project.org
Subject: [R-sig-eco] Comparing multiple distributions within treatments

Dear R Ecologists,
   Does anyone have ideas on how to compare "distributions" (binned count values) between 2 treatments when there are several paired replicates?  I have data on the tolerance to poor water conditions of organisms in two types of streams.  The organisms were sampled from streams that were paired so both treatments were represented in a pair.  I binned the organism counts into bins of tolerance ranges of 1.  Tolerance is on an ordinal by not necessarily ratio scale.  The code below will plot a subset of the data.  What I would like to do is test for a difference between the "distributions" by treatment including the paired set-up.  A Chi-square based test will allow me to compare one distribution of each treatment if I ignore the pairing. Examining the differences within the pairs (e.g., 2stg -- trad for all classes) will compare for one pair.  I think there is a way to combine these approaches to make this comparison, but I am uncertain how to do this.  Perhaps I should get t!
 he differences within the pairs, and then compare the mean "distribution" to a randomized one?  If the treatments are similar, the difference line should be around the diff=0 line, so there may be a test to look at this.
   Any advice appreciated,
   Jeff Holland

my.sites <- c("kirk", "JC", "BG", "vine", "200S", "WS")		# unique site names
my.treat <- c(rep("2st",3),rep("trad",3))					# 2 treatments
my.pair <- c("B", "D", "E", "B", "D", "E")				# 3 sets of paired sites
my.data <- data.frame(row.names=my.sites, t1=c(0,0,0,0,0,0), t2=c(28,16,31,7,0,0),
	t3=c(34,84,5,23,2,12), t4=c(23,68,139,117,2,16) )		# counts within 4 tolerance classes
my.data2 <- cbind(my.treat, my.pair, my.data)
my.data2

attach(my.data2)
plot(seq(0,4,by=4/5), t4, type="none", xlab="Tolerance", ylab="Number") for(d in 1:6){lines(seq(1:4),my.data2[d,3:6], type="b", col=my.data2[d,1], pch=as.numeric(my.data2[d,1]))}
	legend(0.5,120, legend=c("2st","trad"), col=c("black","red"), pch=c(1,2))

_______________________________________________
R-sig-ecology mailing list
R-sig-ecology at r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-ecology



More information about the R-sig-ecology mailing list