|
|
| Line 13: |
Line 13: |
|
| |
|
| <rcode name="answer"> | | <rcode name="answer"> |
| | |
| | library(OpasnetBaseUtils) |
| | library(xtable) |
| | data <- opasnet.csv("5/54/Tigre-data.csv", wiki = "opasnet_en", sep = ",") |
| | data <- data[data$Species != 1 & data$Individuals != 0, ] |
| | head(data) |
| | species <- data$Species |
| | amount <- data$Individuals |
| | transect <- data$Transect |
|
| |
|
| true <- function(pi, q){ | | true <- function(pi, q){ |
| if(q != 1){qD <- (sum(pi^q))^(1 / (1 - q))} else | | qD <- rep(NA, length(q)) |
| {qD <- exp(-sum(pi * log(pi)))} | | for(i in 1:length(q)){ |
| | if(q[i] == 1) |
| | {qD[i] <- exp(-sum(pi * log(pi)))} else |
| | if( q[i] == 999){ |
| | qD[i] <- 1/max(pi)} else |
| | {qD[i] <- (sum(pi^q[i]))^(1 / (1 - q[i]))}} |
| return(qD) | | return(qD) |
| } | | } |
|
| |
|
| | ## JOTAIN KUMMAA TRUE-FUNKTIOSSA, KUN ANTAA JÄRJETTÖMIÄ TULOKSIA |
|
| |
|
| true2 <- function(divj, q){
| | diversity2 <- function(pi){ |
| if(q != 1){qDa <- (sum(wj * divj^(1-q)))^1/(1 - q)} else
| | out <- true(pi, c(0, 1, 2, 999)) |
| {qDa <- exp(sum(wj * log(divj)))}
| | return(out) |
| return(qDa) | |
| } | | } |
|
| |
|
| diversity <- function(amount = rep(1,length(species)), species = 1:length(amount)){
| | diversity2(c(0.2, 0.4, 0.1, 0.1, 0.2)) |
| pi <- tapply(amount, species, sum)
| |
| pi <- pi/sum(pi)
| |
|
| |
|
| out <- data.frame(Name = c("True diversity with q=0", "True diversity with q=1",
| | truelist <- function(pi){ |
| "True diversity with q=2", "True diversity with q=∞", "Richness", "Shannon index",
| | return(true(pi, c(0, 1, 2, 999))) |
| "Simpson index", "Inverse Simpson index", "Gini-Simpson index", "Berger-Parker index"),
| |
| Symbol = c("0D", "1D", "2D", "∞D", "S", "H' or log(1D)", "λ or 1/(2D)", "1/λ or 2D",
| |
| "1-λ or 1-1/(2D)", "1/(∞D)"),
| |
| Value = c(true(pi,0), true(pi,1), true(pi,2), 1/max(pi), length(pi), log(true(pi,1)),
| |
| 1/true(pi,2), true(pi,2), 1-1/true(pi,2), max(pi)))
| |
| return(out)
| |
| } | | } |
|
| |
|
| diversity2 <- function(amount = rep(1,length(species)), species = 1:length(amount)){
| | true2 <- function(divj, wj){ |
| out <- diversity(amount, species)[1:4, 3]
| | #q <- data.frame(Zero = rep(0, length(wj)), One = 1, Two = 2, Infi = 999) |
| return(out) | | if(q == 1) |
| | {qDa <- exp(sum(wj * log(divj)))} |
| | else |
| | {qDa <- (sum(wj * divj^(1-q)))^1/(1 - q)} |
| | return(qDa) |
| } | | } |
|
| |
|
| library(OpasnetBaseUtils)
| | # diversity <- function(amount = rep(1,length(species)), species = 1:length(amount), transect = 1){ |
| library(xtable)
| | |
| data <- opasnet.csv("5/54/Tigre-data.csv", wiki = "opasnet_en", sep = ",")
| | |
| data <- data[data$Species != 1 & data$Individuals != 0, ]
| | # pi <- tapply(amount, species, sum) |
| head(data)
| | # pi <- pi/sum(pi) |
| out <- diversity(species = data$Species, amount = data$Individuals)
| | |
| #print(xtable(out), type = 'html') | |
|
| |
|
| pi <- as.data.frame(as.table(tapply(data$Individuals, data[, c("Transect", "Species")], sum))) | | pij <- as.data.frame(as.table(tapply(amount, data.frame(Transect = transect, Species = species), sum))) |
| pi[is.na(pi)] <- 0 | | pij[is.na(pij)] <- 0 |
| head(pi) | | head(pij) |
| pij <- pi
| |
| colnames(pij)[3] <- "pij" | | colnames(pij)[3] <- "pij" |
|
| |
|
| pij$pij <- pi$Freq/sum(pi$Freq) | | pij$pij <- pij$pij/sum(pij$pij) |
| | pij <- pij[pij$pij != 0, ] |
| | |
| head(pij) | | head(pij) |
| wj <- as.data.frame(as.table(tapply(pij$pij, pij$Transect, sum))) | | wj <- as.data.frame(as.table(tapply(pij$pij, pij$Transect, sum))) |
| Line 66: |
Line 76: |
| colnames(wj) <- c("Transect", "Weight") | | colnames(wj) <- c("Transect", "Weight") |
| out <- merge(pij, wj) | | out <- merge(pij, wj) |
| out <- out[out$pij != 0, ] | | head(out) |
| out <- tapply(out$pij, out$Transect, diversity2) | | out <- tapply(out$pij, out$Transect, truelist) |
| out3 <- data.frame(Zero = 1:nrow(wj), One = NA, Two = NA, Infi = NA)
| | head(out) |
| | qDj <- data.frame(Zero = 1:nrow(wj), One = NA, Two = NA, Infi = NA) |
| for(i in 1:nrow(wj)){ | | for(i in 1:nrow(wj)){ |
| out3[i,] <- out[[i]][1:4]
| | qDj[i,] <- out[[i]][1:4] |
| } | | } |
| q <- c(0, 1, 2, | | qDj |
| out3 <-
| | |
| | #NYT PITÄISI LASKEA true2 MUTTA MITEN YHDISTÄNKÄÄN div ja q? |
| | |
| | out <- data.frame(Name = c("True diversity with q=0", "True diversity with q=1", |
| | "True diversity with q=2", "True diversity with q=∞", "Richness", "Shannon index", |
| | "Simpson index", "Inverse Simpson index", "Gini-Simpson index", "Berger-Parker index"), |
| | Symbol = c("0D", "1D", "2D", "∞D", "S", "H' or log(1D)", "λ or 1/(2D)", "1/λ or 2D", |
| | "1-λ or 1-1/(2D)", "1/(∞D)"), |
| | Value = c(true(pi,0), true(pi,1), true(pi,2), 1/max(pi), length(pi), log(true(pi,1)), |
| | 1/true(pi,2), true(pi,2), 1-1/true(pi,2), max(pi))) |
| | #return(out) |
| | #} |
| | |
| | print(xtable(out), type = 'html') |
|
| |
|
|
| |
|
Question
How to calculate diversity indices?
Answer
Upload your data to Opasnet Base. Use the function diversity to calculate the most common indices.
Actual function diversity
+ Show code- Hide code
library(OpasnetBaseUtils)
library(xtable)
data <- opasnet.csv("5/54/Tigre-data.csv", wiki = "opasnet_en", sep = ",")
data <- data[data$Species != 1 & data$Individuals != 0, ]
head(data)
species <- data$Species
amount <- data$Individuals
transect <- data$Transect
true <- function(pi, q){
qD <- rep(NA, length(q))
for(i in 1:length(q)){
if(q[i] == 1)
{qD[i] <- exp(-sum(pi * log(pi)))} else
if( q[i] == 999){
qD[i] <- 1/max(pi)} else
{qD[i] <- (sum(pi^q[i]))^(1 / (1 - q[i]))}}
return(qD)
}
## JOTAIN KUMMAA TRUE-FUNKTIOSSA, KUN ANTAA JÄRJETTÖMIÄ TULOKSIA
diversity2 <- function(pi){
out <- true(pi, c(0, 1, 2, 999))
return(out)
}
diversity2(c(0.2, 0.4, 0.1, 0.1, 0.2))
truelist <- function(pi){
return(true(pi, c(0, 1, 2, 999)))
}
true2 <- function(divj, wj){
#q <- data.frame(Zero = rep(0, length(wj)), One = 1, Two = 2, Infi = 999)
if(q == 1)
{qDa <- exp(sum(wj * log(divj)))}
else
{qDa <- (sum(wj * divj^(1-q)))^1/(1 - q)}
return(qDa)
}
# diversity <- function(amount = rep(1,length(species)), species = 1:length(amount), transect = 1){
# pi <- tapply(amount, species, sum)
# pi <- pi/sum(pi)
pij <- as.data.frame(as.table(tapply(amount, data.frame(Transect = transect, Species = species), sum)))
pij[is.na(pij)] <- 0
head(pij)
colnames(pij)[3] <- "pij"
pij$pij <- pij$pij/sum(pij$pij)
pij <- pij[pij$pij != 0, ]
head(pij)
wj <- as.data.frame(as.table(tapply(pij$pij, pij$Transect, sum)))
head(wj)
colnames(wj) <- c("Transect", "Weight")
out <- merge(pij, wj)
head(out)
out <- tapply(out$pij, out$Transect, truelist)
head(out)
qDj <- data.frame(Zero = 1:nrow(wj), One = NA, Two = NA, Infi = NA)
for(i in 1:nrow(wj)){
qDj[i,] <- out[[i]][1:4]
}
qDj
#NYT PITÄISI LASKEA true2 MUTTA MITEN YHDISTÄNKÄÄN div ja q?
out <- data.frame(Name = c("True diversity with q=0", "True diversity with q=1",
"True diversity with q=2", "True diversity with q=∞", "Richness", "Shannon index",
"Simpson index", "Inverse Simpson index", "Gini-Simpson index", "Berger-Parker index"),
Symbol = c("0D", "1D", "2D", "∞D", "S", "H' or log(1D)", "λ or 1/(2D)", "1/λ or 2D",
"1-λ or 1-1/(2D)", "1/(∞D)"),
Value = c(true(pi,0), true(pi,1), true(pi,2), 1/max(pi), length(pi), log(true(pi,1)),
1/true(pi,2), true(pi,2), 1-1/true(pi,2), max(pi)))
#return(out)
#}
print(xtable(out), type = 'html')
#pi.j <-
out <- data.frame(Name = c("True diversity with q=0", "True diversity with q=1",
"True diversity with q=2", "True diversity with q=∞", "Richness", "Shannon index",
"Simpson index", "Inverse Simpson index", "Gini-Simpson index", "Berger-Parker index"),
Symbol = c("0D", "1D", "2D", "∞D", "S", "H' or log(1D)", "λ or 1/(2D)", "1/λ or 2D",
"1-λ or 1-1/(2D)", "1/(∞D)"),
Value = c(true(pi,0), true(pi,1), true(pi,2), 1/max(pi), length(pi), log(true(pi,1)),
1/true(pi,2), true(pi,2), 1-1/true(pi,2), max(pi)))
out
| |
Example to use function
+ Show code- Hide code
library(xtable)
if(is.null(data)) {data <- c(1,4,6,9,3,4,5,6,5,4,3,3,5,5,7,5,5,4,3,4,5,6,8,9,5,4,5,4,3,3,4,9,6,6,4,5,3,2,1,1,2,3,4,3,2,3,4,5,6,7,7)}
data
if(individual==TRUE) out <- diversity(species = data) else out <- diversity(amount = data)
print(xtable(out), type = 'html')
| |
The data should be given in R format as a list of values in parenthesis, beginning with c:
c(3,5,3,5,2,1,3,3,4,2) or equivalently c(0.1,0.2,0.4,0.1,0.2)
where the values are either
- identifiers of the species 1,2,3... in which the individuals belong (one entry per individual), or
- abundancies of species, i.e. proportions of individuals belonging to each species among the whole population (one entry per species).
Rationale
Diversity indices are thoroughly described in Wikipedia.
See also
References
Related files
<mfanonymousfilelist></mfanonymousfilelist>