Os métodos de classificação agrupam as parcelas conforme a similaridade na composição em grupos. Isso pode ser muito interessante para verificarmos se há conjuntos discretos de objetos, no nosso caso parcelas. Para tanto precisamos primeiro calcular um índice que nos diga o quanto cada parcela é similar as outras. A primeira parte desse exercicio é como calcular medidas de distância e similaridade na nossa amostra em seguida usamos um método de agrupamento dessas parcelas para finalmente apresentar os dados graficamente em um dendrograma.

Existem uma infinidade de índices utilizados na ecologia para medir a similaridade ou dissimilaridade (distância) entre objetos. Para a comparação entre parcelas podemos ter normalmente dois tipos de dados referentes às espécies: presença ou abundância. Quando a natureza da informação é apenas a presença fazemos uso de indices quantitativos, o mais intuitivo é o de Jaccard, que nos diz a proporção de espécies que é compartilhada pelas parcelas.

Para calcular a similaridade de Jaccard vamo usar a formula:

$$Sim_(jacc)= (nsp_(com))/(nsp_(1) + nsp_(2)- nsp_(com))

Vamos usar as seguintes parcelas hipotéticas para calcular os índices:

Siga os passos abaixo para o cálculo no R!

 
plot1=c(2,13,0,1, 20)
plot1
plot2=c(8,6,5,15,0)
plot2
amostra=data.frame(plot1,plot2)
nsp_1=sum(plot1>0)
nsp_1
nsp_2=sum(plot2>0)
nsp_2
nsp_com=sum(plot1>0 & plot2>0)
nsp_com
jacc= nsp_com/(nsp_1 + nsp_2 - nsp_com)
jacc
## criando o objeto amostra com as duas parcelas:
amostra=data.frame(plot1,plot2)

Temos portanto, 60% de similaridade entre a parcela 1 e 2 segundo o índice de Jaccard.

Quando queremos comparar não só a presença qualitativa das espécies, mas suas abundâncias precisamos usar um índice quantitativo. Uma das maneiras de pensar nesses indices é pela medida de distância, ou seja pelo quanto uma parcela se diferencia da outra medida pela diferença entre as abundâncias das espécies. Caso esteja interessado veja a explicação da distância euclidiana que é bastante intuitiva geometricamente. No nosso caso vamos usar um índice conhecido como Bray-Curtis. sendo:

no. de indivíduos da espécie i na parcela1: $$nsp_i_1

no. de indivíduos da espécie i na parcela: $$nsp_i_2

S= número de espécies total Temos a distância de Bray Curtis dada por:

$$BC_(dist)=(sum_(n=1)^S(|nsp_(i1)- nsp_(i2)|))/(sum_(n=1)^S(nsp_(i1)+ nsp_(i2)))

e sua similaridade por:

$$BC_(sim)=(sum_(n=1)^S(min(nsp_(i1), nsp_(i2))))/(sum_(n=1)^S(nsp_(i1)+ nsp_(i2)))

Vamos calcular então similaridade para as duas parcela:

bc.sim=sum(apply(amostra, 1, min))/(sum (amostra$plot1+amostra$plot2))
bc.sim

Por que os valores dos índices são tão diferentes?

Agora vamos aplicar esse mesmo índice para nossas comunidades discretas e contínuas. Para tanto temos que fazer uma função que compare cada uma das parcelas com as outras todas.Ei-la:

sim<-function(dados, indice="bc")
	{
	nplot=dim(dados)[2]
	similar=matrix(1,ncol=nplot,nrow=nplot)
	rownames(similar)<-paste("plot", c(1:nplot))
	colnames(similar)<-paste("plot", c(1:nplot))
		for(i in 1:(nplot-1))
		{
		m=i+1
		for(m in m:nplot)
		{
		if(indice=="jacc")
			{
			dados[dados>0]=1
			co.oc=sum(dados[,i]>0 & dados[,m]>0)
			total.sp=sum(dados[,i])+sum(dados[,m])-co.oc
			similar[i,m]=co.oc/total.sp 
			similar[m,i]=co.oc/total.sp
			}
		if(indice=="bc") 
			{
			bc.sim=sum(apply(dados[,c(i,m)], 1, min))/(sum (dados[,c(i,m)]))
			similar[i,m]=bc.sim
			similar[m,i]=bc.sim
			}
		}
		}
	return(round(similar,3))
	}

Agora vamos aplicar essa função para a nossa comunidade virtual contínua:

sim.cont1=sim(com1.cont, indice="bc")
sim.cont1
round(sim.cont1, 3)
sim.cont1[1,4]
sim.cont1[4,1]
sim.cont1[5,7]
sim.cont1[7,5]

O que significam os valores 1 nesse resultado? Por que os valores na posição linha/coluna é igual ao coluna/linha? Consegue identificar padrões nesses dados? Faça o mesmo para a comunidade discreta, mas guarde o resultado em um objeto de nome sim.disc1. Chegamos nas matrizes de similiaridades, o primeiro passo da nossa classificação.

Há vários métodos e algorítmos¹⁾ para agrupar os objetos em uma análise de classificação. Primeiro precisamos decidir se queremos iniciar juntando os elementos (aglomerativa) ou se queremos tomar todo o conjunto de objetos e ir separando em grupos (divisiva). No nosso caso vamos iniciar com as parcela como unidades e vamos montando grupos sucessivamente a partir daqueles já formados. Isso se chama análise de classificação (ou agrupamento) hierarquica aglomerativa. Apesar do nome feio, não há nada de complicado na lógica da análise. Vamos acompanhar passo a passo.

• 1. Vamos procurar entre os pares de parcelas, aquela que tem mais similaridade

max(sim.cont1)
sim.cont1==max(sim.cont1)

O verdadeiro (TRUE) significa a posição onde os valores máximos se encontram. Essas similaridades não nos interessa. Não há porque agrupar algo com ele mesmo. Portanto, vamos retira do nossa matriz de similaridade essa informação assim com a redundância.

nome.par=rownames(sim.cont1)
upper.tri(sim.cont1, diag=TRUE)
sim.cont1[upper.tri(sim.cont1, diag=TRUE)]=NA
sim.cont1

Vamos agora procurar o maior valor e depois saber qual o par de parcelas que apresenta essa similiaridade:

max1=max(sim.cont1,na.rm=TRUE)
max1
maior1=which(sim.cont1==max(sim.cont1,na.rm=TRUE), arr.ind=TRUE)
maior1

Vamo tomar apena o primeiro valor maior, caso tenham outros.

par1=nome.par[maior1[1,]]
par1
cat("\n\t 1a. ligação ", paste(par1[1],par1[2], sep=" x "), "; ligação = ", max1, "\n" )

Esse é nosso primeiro grupo: Agora temos que tomar uma outra decisão: como esse grupo agora formado irá se ligar aos outras parcelas? Podemos decidir que ele irá se ligar aos outros pelo valor da máxima similaridade de seus componentes (ligação máxima), ou pela minima. Vamos usar a média do grupo:

g1.n1=strsplit(nome.par[maior1[1,2]]," ")[[1]][2]
g1.n2=strsplit(nome.par[maior1[1,1]]," ")[[1]][2]
g1.nome=paste("g", paste(g1.n1,g1.n2, sep=","))
g1.nome

mat.g1=sim.cont1[-maior1[1,],-maior1[1,]]
g1a=apply(sim.cont1[maior1[1,],-maior1[1,]],2,mean)
g1a[is.na(g1a)]=0
g1b=apply(sim.cont1[-maior1[1,],maior1[1,]],1,mean)
g1b[is.na(g1b)]=0
gr1=rbind(mat.g1,g1a+g1b)
grupo1=cbind(gr1,NA)
rownames(grupo1)<-c(nome.par[-maior1[1,]],g1.nome)
colnames(grupo1)[dim(grupo1)[2]]<-g1.nome
grupo1

Agora repetimos os mesmos passos procurando o novo grupo,nessa nova matriz:

nome.par2=rownames(grupo1)
max2=max(grupo1,na.rm=TRUE)
max2
maior2=which(grupo1==max(grupo1,na.rm=TRUE), arr.ind=TRUE)
maior2
g2.n1=strsplit(rownames(grupo1)[maior2[1,2]]," ")[[1]][2]
g2.n2=strsplit(rownames(grupo1)[maior2[1,1]]," ")[[1]][2]
g2.nome=paste(paste("g",g2.n1,sep="_"),g2.n2, sep=",")  
g2.nome
cat("\n\n\t 2a. ligação ", paste(nome.par2[maior2[1,2]],nome.par2[maior2[1,1]], sep=" x "), "; ligação = ", max2, "\n" )

mat.g2=grupo1[-maior2[1,],-maior2[1,]]
g2a=apply(grupo1[maior2[1,],-maior2[1,]],2,mean)
g2a[is.na(g2a)]=0
g2b=apply(grupo1[-maior2[1,],maior2[1,]],1,mean)
g2b[is.na(g2b)]=0
gr2=rbind(mat.g2,g2a+g2b)
grupo2=cbind(gr2,NA)
rownames(grupo2)<-c(nome.par2[-maior2[1,]],g2.nome)
colnames(grupo2)[dim(grupo2)[2]]<-g2.nome
grupo2

Agora mais uma vez buscamos a maior ligação entre todas na nova matriz e formamos novo grupo.

nome.par3=rownames(grupo2)
max3=max(grupo2,na.rm=TRUE)
max3
maior3=which(grupo2==max(grupo2,na.rm=TRUE), arr.ind=TRUE)
maior3
g3.n1=strsplit(rownames(grupo2)[maior3[1,2]]," ")[[1]][2]
g3.n2=strsplit(rownames(grupo2)[maior3[1,1]]," ")[[1]][2]
g3.nome=paste(paste("g",g3.n1,sep="_"),g3.n2, sep=",")  
g3.nome
cat("\n\n\t 3a. ligação ", paste(nome.par3[maior3[1,2]],nome.par3[maior3[1,1]], sep=" x "), "; ligação = ", max3, "\n" )

mat.g3=grupo2[-maior3[1,],-maior3[1,]]
g3a=apply(grupo2[maior3[1,],-maior3[1,]],2,mean)
g3a[is.na(g3a)]=0
g3b=apply(grupo2[-maior3[1,],maior3[1,]],1,mean)
g3b[is.na(g3b)]=0
gr3=rbind(mat.g3,g3a+g3b)
grupo3=cbind(gr3,NA)
rownames(grupo3)<-c(nome.par3[-maior3[1,]],g3.nome)
colnames(grupo3)[dim(grupo3)[2]]<-g3.nome
grupo3

Se continuarmos a fazer a mesma operação teremos um único grupo no final onde existe a última ligação entre os grupos formados anteriormente. Parece bastante tedioso! Por sorte temos uma função no pacote básico do R que da conta do recado para nós. Entendendo o que está sendo feito, podemos usar a função. Caso tenha ainda alguma dúvida sobre o processo de classificação hierarquica aglomerativa, refaça os passos anteriores ou peça socorro!. Caso tenha entendido já pode perguntar numa rodinha na Vila Madalena: - “Vocês conhecem o algorítmo hiearquico aglomerativo?”. Seu sucesso reprodutivo vai ficar em baixa, mas vale a pena para ver a cara dos coleguinhas… Vamos agora rodar com a nossa função de similaridade e com a função básica do R hclust():

  clas.cont1=hclust(as.dist(1-sim.cont1), method="average")
  dend.cont1=as.dendrogram(clas.cont1, hang=-1)
  plot(dend.cont1)

Vamos testar comparar os outros métodos de ligação:

par(mfrow=c(2,2))
clas.cont1a=hclust(as.dist(1-sim.cont1), method="single")
plot(as.dendrogram(clas.cont1, hang=-1), main="Ligação simples",ylab="Bray-Curtis")
clas.cont1b=hclust(as.dist(1-sim.cont1), method="complete")
plot(as.dendrogram(clas.cont1b, hang=-1), ylab="Bray-Curtis",main="Ligação completa") 
clas.cont1c=hclust(as.dist(1-sim.cont1), method="average")
plot(as.dendrogram(clas.cont1c, hang=-1), ylab="Bray-Curtis",main="Ligação média")
clas.cont1d=hclust(as.dist(1-sim.cont1), method="centroid")
plot(as.dendrogram(clas.cont1d, hang=-1),ylab="Bray-Curtis",main="Ligação centroide")

Tente interpretar o que está acontecendo nas diferentes ligações. Para saber mais sobre os métodos de ligação consulte http://en.wikipedia.org/wiki/Hierarchical_clustering e http://en.wikipedia.org/wiki/Centroid

Agora é por sua conta. Vamos usar as suas comunidades virtuais para ver como elas são classificadas. Abaixo apresentamos o código básico a seguir. Onde está escrito “minha.com” deve colocar o nome do objeto onde guardou a amostra da sua comunidade virtual contínua e discreta. Troque também os nomes que tenha “*. minha” para guardar os resultados, note que esses nomes de objetos não devem ter “ ”, enquanto alguns argumentos necessariamente precisam.

sim.minha=sim(**minha.com**, indice="bc")
clas.minha=hclust(as.dist(1-sim.minha), method="average")
plot(as.dendrogram(clas.minha, hang=-1), ylab="Bray-Curtis",main="Ligação média")

Salve os gráficos das comunidades e interprete as diferenças entre as suas comunidades: continua e discreta. Para fazer gráficos em um mesmo dispositivo use o comando:

par(mfrow=c(2,2))

O primeiro valor significa o numero de linhas e o segundo o número de colunas em que sua janela gráfica será dividida. Para retornar ao padrão inicial de apenas um gráfico por janela, feche o dispositivo para iniciar outro no momento que a função de gráfico for chamada.