Cantinho do R

Curvas de acúmulo de espécies bonitinhas

Não gosto daquele gráfico com cara de banana do vegan!!!

Por Marcos Vital, do LEQ-UFAL

Este material foi construído com a ajuda de muitas pessoas que acreditam no LEQ e em Ciência Livre. Muito obrigado!

Para mais material, visite o Cantinho do R


Recentemente, durante uma disciplina aqui da UFAL (a Feramentas Computacionais em Ecologia), eu e os alunos estávamos estudando as curvas de acúmulo de espécie, e brincando com as funções do pacote vegan e o método de rarefação. Até aí tudo bem, tudo muito divertido (juro, é bem legal!), mas rolou um problema… :P

É que o gráfico padrão do vegan pra essas curvas é, tadinho, super feioso. Olha ele aqui, ó:

library(vegan)
## Loading required package: permute
## Loading required package: lattice
## This is vegan 2.4-1
data(dune)

resultado<-specaccum(dune)

plot(resultado)

Bleh!!!!! Não gostei…

Daí começamos a escrever um script bacaninha, que pudesse fazer uma curva legalzinha. A coisa começou a andar, o novo gráfico ficou super bonitinho, mas nesta disciplina o objetivo era puxar um pouco as coisas, então propus: vamos transformar nosso script em uma função, de forma que a gente possa automatizar a construção do gráfico! Bora? Bora!!!!

E fizemos. Agora, neste tutorialzinho, quero apresentar e disponibilizar a função para vocês usarem, curtirem, modificarem e aproveitarem como acharem melhor. ;)

 


 

Criando a função - começando entendendo o resultado gerado no vegan

Vamos lá. Primeiro, o nosso objetivo: criar uma função que leia o resultado da função specaccum do vegan. Claro que você tem que sacar o que a função faz, o que não é o objetivo deste tutorial, beleza? Em outro momento posso explorar mais coisas sobre isso, mas neste momento pressuponho que você sabe usar a função.

Se vamos usar a função, temos que entender o que tem “dentro” do objeto que ela gera. O truque aqui, super útil em outras situações similares, é usar a função str para olhar o objeto por dentro, assim, ó:

str(resultado)
## List of 7
##  $ call    : language specaccum(comm = dune)
##  $ method  : chr "exact"
##  $ sites   : int [1:20] 1 2 3 4 5 6 7 8 9 10 ...
##  $ richness: num [1:20] 9.85 15.11 18.51 20.94 22.75 ...
##  $ sd      : num [1:20] 2.35 1.88 1.57 1.45 1.39 ...
##  $ perm    : NULL
##  $ freq    : Named num [1:30] 7 10 2 8 6 6 5 1 1 2 ...
##   ..- attr(*, "names")= chr [1:30] "Achimill" "Agrostol" "Airaprae" "Alopgeni" ...
##  - attr(*, "class")= chr "specaccum"

Olha que lindeza, né não?!

O que o R está nos mostrando são as “partes” do objeto que chamamos de resultado, e que foi gerado pela specaccum. Os pedacinhos que nos interessão são os valores de riqueza e o erro padrão: richness e sd (sim, sd é o erro padrão neste caso; looooonga história, fica pra outra hora). Ah sim, e podemos pegar a medida de esforço amostral, guardada em samples, pois servirá como nosso eixo horizontal.

Quer conferir? Fuce um pouco, é sempre saudável fuçar as coisas assim no R:

resultado$richness
##  [1]  9.85000 15.11053 18.51053 20.93746 22.75432 24.14959 25.23956
##  [8] 26.10355 26.79824 27.36496 27.83398 28.22755 28.56204 28.84964
## [15] 29.09959 29.31909 29.51404 29.68947 29.85000 30.00000
resultado$sd
##  [1] 2.3510636 1.8763851 1.5722711 1.4469584 1.3901594 1.3530349 1.3164796
##  [8] 1.2749034 1.2282010 1.1763410 1.1193437 1.0564537 0.9874094 0.9115998
## [15] 0.8286890 0.7380921 0.6333903 0.5139710 0.3570714 0.0000000
resultado$Samples
## NULL

Se quiser, confira o que tem dentro de tudo, é um aprendizado legal. Se você não sabia, agora sabe que a função str é uma super mão na roda.

Agora conhecemos os pedacinhos, e podemos seguir adiante.

 

Criando a função - começando entendendo o resultado gerado no vegan

Bom, podemos começar salvando cada pedacinho num objeto. Isso não é verdadeiramente necessário (a gente pode continuar usando resultado$pedaço-que-queremos), mas acredito que fique um pouco mais fácil de acompanhar a coisa toda, e os nomes ficam menores. Então vamos lá:

#Pegando os resultados que queremos:

riqueza<-resultado$richness
desvio<-resultado$sd
amostras<-resultado$sites

Beleza, cada elemento é um objeto com nome fácil de lembrar e escrever. ;)

 

Criando a função - preparação da escala do gráfico

Então, o nosso gráfico vai começar com a curva dos valores de riqueza. Até aí tudo bem, mas tem um probleminha: como além da riqueza vamos incluir o erro padrão, o R pode acabar criando uma escala no eixo y que não alcança os valores necessários para as barras de erro, e pode acabar “comendo” um pedaço do gráfico.

Resulver isso é simples, basta a gente estabelecer os limites do gráfico com o bom e velho argumento ylim. Mas aí tem uma coisa: eu quero que isso seja feito automaticamente, pra não ficarmos tendo que estabelecer isso manualmente toda santa vez…

Então vamos usar os próprios valores de riqueza como guia, e depois “espichar” os limites pra um pouco mais que a maior riqueza e um pouco menos que a menor riqueza. Sacou?

Fazendo beeeeem passo a passo, pra ficar fácil de entender, o esquema é:

#Achando o valor máximo e mínimo:
y.mínimo<-min(riqueza)
y.máximo<-max(riqueza)

#Estabelecendo valores um pouco abaixo do mínimo e um pouco acima do máximo:
y.baixo<-y.mínimo*0.80
y.alto<-y.máximo+(y.mínimo*0.20)

#Salvando tudo num objeto que vai guardar estes limites:
limites.y<-c(y.baixo, y.alto)

Ufa, pronto. Agora teremos um objeto, o limites.y que guarda a informação do tamanho do eixo y. Prático, né?

 

Criando a função - agora sim, o gráficoooooooooo

Bora, bora, bora!

Vamos fazer o gráfico em três etapas: curva, barras de erro e pontos. Lá vai:

#Criando a curva:
plot(riqueza~amostras, type="l", ylim=limites.y, las=1, xlab="Esforço amostral", ylab="Riqueza de espécies")

#Adicionando as barras de erro:
arrows(amostras, riqueza-desvio, amostras, riqueza+desvio, angle=90, code=3, length=0.05)
## Warning in arrows(amostras, riqueza - desvio, amostras, riqueza + desvio, :
## zero-length arrow is of indeterminate angle and so skipped
#E adicionando os pontos:
points(riqueza~amostras, pch=21, bg="white")

Nóóóóó, ficou bonitinho, não foi não? :D (o nóóó foi o meu espírito mineiro falando)

Quer saber mais sobre a função arrows? Veja este script aqui, nele eu explico.

Mas falta uma coisaaaaa!!!!11!um!!!onze!!!!!

hehe

Falta o seguinte: bora colocar essa paradinha toda dentro de uma função! Assim fica tudo mais fácil, a gente “chama” a função e pronto. Aqui como a função ficou, completinha:

curvacum<-function(resultado, nome.x, nome.y) {
    riqueza<-resultado$richness
    desvio<-resultado$sd
    amostras<-resultado$sites

    y.mínimo<-min(riqueza)
    y.máximo<-max(riqueza)
    y.baixo<-y.mínimo*0.80
    y.alto<-y.máximo+(y.mínimo*0.20)
    limites.y<-c(y.baixo, y.alto)

    plot(riqueza~amostras, type="l", ylim=limites.y, las=1, xlab=nome.x, ylab=nome.y)
    arrows(amostras, riqueza-desvio, amostras, riqueza+desvio, angle=90, code=3, length=0.05)
    points(riqueza~amostras, pch=21, bg="white")

}

Do jeitinho que está, a nossa função curvacum pede o seguinte: o resultado da função specaccum do pacote vegan, o nome do eixo x e o nome do eixo y. Legal, né?

Se você executar o código acima, a função agora está na memória do R, e aí você pode fazer assim (com outros dados agora):

data(BCI)

resultado2<-specaccum(BCI)

curvacum(resultado2, "Esforço amostral", "Riqueza de espécies")

Maneiríssimo, né? Mas pode ficar ainda mais prático!

Se você salvar o script com a função, pode “chamá-lo” a qualquer hora com a função source. DAí basta você copiar o script, colar na sua pasta de trabalho do momento, usar a função e ser feliz (sem ter que incluir toda a função cada vez que usar). Bom, né?

Então você pode baixar o script prontinho aqui, e depois fazer assim, ó:

#Carrega a função pronta:
source("curvacum.R")

#Vamos usar ainda outro exemplo?
data(mite)

#Vegan, eu escolho você!
resultado3<-specaccum(mite)

#curvacum, vai!
curvacum(resultado3, "Esforço amostral", "Riqueza de espécies")

E é isso aí. Espero que vocês tenham se divertido. Ou pelo menos gstado um pouquinho. Ou pelo menos usem a função depois, sei lá.

Um abração, e até a próxima!

Prof Marcos

 


Este documento faz parte do material que disponibilizo no meu blog, o “Cantinho do R”, e foi foi criado utilizando o software R, o programa RStudio e a linguagem Markdown. Para saber mais sobre eles, acesse:

http://www.r-project.org/

http://www.rstudio.com

http://rmarkdown.rstudio.com/

 

Quer comentar ou perguntar? Faça isso no blog, aqui