[R] Some code to run Hutcheson t-test using R
Luigi Marongiu
m@rong|u@|u|g| @end|ng |rom gm@||@com
Thu Sep 10 10:35:45 CEST 2020
Thank you very much for the code, that was very helpful.
I got the article by Hutcheson -- I don't know if I can distribute it
, given the possible copyrights, or if I can attach it here -- but it
does not report numbers directly: it refers to a previous article
counting bird death on a telegraph each year. The numbers
are:
bird_1956 <- c(4,4,190,135,56,3,2,2,1,12,41,201,1,0,131,3)
bird_1957 <- c(4,111,53,66,146,222,126,61,0,2323,21)
bird_1958 <- c(0,3,32,228,56,102,0,11,2,220,0)
bird_1959 <- c(0,0,14,59,26,68,0)
bird_1960 <- c(0,0,73,66,71,25,0,109,63,1)
This for sake of the argument.
As for my problem, I implemented the Shannon index with the package
iNext, which only gives me the index itself and the 95% CI. Even when
I implemented it with vegan, I only got the index. Essentially I don't
have a count of species I could feed into the Hutcheson's. Is there a
way to extract these data? Or to run a Hutcheson's on the final index?
Thank you
On Tue, Sep 8, 2020 at 7:43 PM Karl Schilling
<karl.schilling using uni-bonn.de> wrote:
>
> Dear Luigi,
>
> below some code I cobbled together based on the Hutcheson paper you
> mentioned. I was lucky to find code to calculate h and, importantly, its
> variance in the R-package QSutils - you may find it on the Bioconductor
> website.
>
> here is the code, along with an example. I also attach the code as an
> R-file.
>
> Hope that helps.
> All my best
>
> Karl
> PS don't forget to adjust for multiple testing if you compare more than
> two groups.
> K
>
>
> # load package needed
> # QSutils is on Bioconductor
> library(QSutils)
>
> # here some exemplary data - these are the data from Pilou 1966 that are
> used
> # in the second example of Hutcheson, J theor Biol 129:151-154 (1970)
>
> earlier <- c(0,0,146,0,5,46,0,1,295,0,0,3,0,0,0,0,0)
> later <- c(0,0,142,0,5,46,0,1,246,0,0,3,0,0,0,0,0)
> # numbers of the first example used by Hutcheson were unfortunately not
> # available to me
>
> # here starts the code ; you may replace the variables "earlier" and "later"
> # by your own numbers.
>
> # calculate h, var(h) etc
> h1 <- Shannon(earlier)
> varh1 <- ShannonVar(earlier)
> n1 <- sum (earlier)
> h2 <- Shannon(later)
> varh2 <- ShannonVar(later)
> n2 <- sum (later)
> degfree <- (varh1 + varh2)**2 /(varh1**2/n1 + varh2**2 /n2)
>
> # compare numbers with those in the paper
> h1
> h2
> varh1
> varh2
> # I assume that minor numerical differences are due to differences in the
> # numerical precision of computers in the early seventies and today / KS
>
> # this is the actual t-test
> t <- (h1-h2) /sqrt(varh1 + varh2)
> p <- 2*pt(-abs(t),df= degfree)
> p
>
> # that's it
> # Best
> # Karl
> --
> Karl Schilling, MD
> Professor of Anatomy and Cell Biology
> Anatomisches Institut
> Rheinische Friedrich-Wilhelms-Universität
> Nussallee 10
>
> D-53115 Bonn
> Germany
>
> phone ++49-228-73-2602
>
--
Best regards,
Luigi
More information about the R-help
mailing list