NOTA: Tem que carregar o “workspace” (arquivo .RData) do exercicio anterior pois lá tem objetos que serão usados nesse exercicio.

Carregar os seguintes pacotes:

library(raster)
library(maptools)
library(rgeos)
library(sp)

Nesse exemplo, a gente vai usar os dados dos morcegos da IUCN disponíveis na internet (http://www.iucnredlist.org/technical-documents/spatial-data) Vcs podem usar qualquer outro grupo dos disponíveis na IUCN ou em outro site (no caso, tem que baixar os “shapefiles”, que são arquivos com extensão ‘.shp’ e os associados)

bats <- readShapePoly("phyllostomidae_shapes/phylllostomidae_marcelo.shp")
plot(bats,col=rainbow(46, alpha=0.5))

“Projetar” em coordenadas geográficas – definir a projeção dos dados/shapefiles

proj4string(bats) <- CRS("+proj=longlat +datum=WGS84")

Gerar um raster para America

amer_ras <- raster()
# estabelecer o "extent" do raster pra ficar justamente nas coordenadas da America
extent(amer_ras) <- c(-150,-20,-60,40)

Se olham o “objeto” gerado, ele teria que ser assim:

amer_ras
class       : RasterLayer 
dimensions  : 180, 360, 64800  (nrow, ncol, ncell)
resolution  : 0.3611111, 0.5555556  (x, y)
extent      : -150, -20, -60, 40  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0

Mudar a resolução do raster da America pra ficar em 1 grau long-lat.

res(amer_ras) <- 1

Agora teria que ser assim

amer_ras
class       : RasterLayer 
dimensions  : 100, 130, 13000  (nrow, ncol, ncell)
resolution  : 1, 1  (x, y)
extent      : -150, -20, -60, 40  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0

Raster de riqueza de espécies

Gerar um raster pra os dados dos morcegos baseado no raster da America já pronto, esse aqui sería o mapa de riqueza de especies de morcegos em America!

Nota: Tem outras formas de gerar o raster de riqueza que são, mais úteis para fazer análises biogeográficos/macroecológicos (vejam, por exemplo, http://bvilela.weebly.com/how-to-make-a-presence-absence-matrix.html)

bats_raster <- rasterize(bats, amer_ras, field="binomial",fun=function(x,...){length(unique(na.omit(x)))})

Graficar o raster dos morcegos.

plot(bats_raster)
plot(wrld_simpl,add=T)

E aí, o mapa ficou legal?

Dependência da escala

Gerar outro raster de America com uma resolução diferente, de 2 graus long-lat

amer_ras2 <- amer_ras
# definir a resolução de 2 graus
res(amer_ras2) <- 2
# o "objeto" tem que ser assim:
amer_ras2
class       : RasterLayer 
dimensions  : 50, 65, 3250  (nrow, ncol, ncell)
resolution  : 2, 2  (x, y)
extent      : -150, -20, -60, 40  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0

Gerar o raster dos morcegos pra nova resolução (usando o novo raster da America em 2 graus long-lat)

bats_raster2 <- rasterize(bats, amer_ras2, field="binomial",fun=function(x,...){length(unique(na.omit(x)))})

Graficar o novo raster.

plot(bats_raster2)
plot(wrld_simpl,add=T)

E agora, ficou bom? é igual ao anterior?

Nota: Pra manter as diferentes graficas dos mapas com diferentes resoluções, tem que abrir uma nova janela de graficação no R. Olhem no manual do R como fazer isso. De outra forma, o grafico será gerado na mesma janela e vão “perder” o anterior. Outra maneira, seria salvar (save as…) o grafico e depois abrir todos eles em algum outro programa pra comparar-los visualmente.

Outra resolução, mais grossa? gerar um raster America com resolução de 4 graus long-lat

amer_ras4 <- amer_ras2

res(amer_ras4) <- 4

Gerar o raster dos morcegos com o novo de America 4 graus long-lat

bats_raster4 <- rasterize(bats, amer_ras4, field="binomial",fun=function(x,...){length(unique(na.omit(x)))})

Graficar.

plot(bats_raster4)
plot(wrld_simpl,add=T)

O padrão mudou?

Resolução mais grossa ainda? gerar novo raster da America em 8 graus long-lat

amer_ras8 <- amer_ras4

res(amer_ras8) <- 8

E o raster de riqueza de morcegos na nova resolução

bats_raster8 <- rasterize(bats, amer_ras8, field="binomial",fun=function(x,...){length(unique(na.omit(x)))})

Graficar.

plot(bats_raster8)
plot(wrld_simpl,add=T)

E aí, mudou? aonde fica mais diferente?

par(mfrow=c(2,2))

plot(bats_raster)
plot(wrld_simpl,add=T)

plot(bats_raster2)
plot(wrld_simpl,add=T)

plot(bats_raster4)
plot(wrld_simpl,add=T)

plot(bats_raster8)
plot(wrld_simpl,add=T)

dev.off()

Agora pro Brasil

Gerar um raster novo (vazio)

ras_2 <- raster()

#definir valores de zero pra o raster
values(ras_2) <- 0
# definir o "extent" pra corresponder com as coordenadas do Brasil
extent(ras_2) <- c(-75,-35,-40,5)
# definir a resolução (por enquanto 1 grau long-lat)
res(ras_2) <- 1

Gerar o raster do Brasil usando o raster gerado antes cortado com o mapa do Brasil que geramos nas anteriormente (“brazil_map”)

ras_brasil <- rasterize(brazil_map,ras_2)

Graficar o novo raster pra o Brasil

plot(ras_brasil)
plot(brazil_map,add=T)

Ficou legal?

Gerar o raster de riqueza de morcegos no Brasil usando o raster do país

bats_brasraster <- rasterize(bats, ras_brasil, field="binomial",fun=function(x,...){length(unique(na.omit(x)))})

Graficar o resultado. Como está distribuida a riqueza de morcegos no Brasil? onde tem mais/menos espécies?

#Gerar uma "máscara" para que o raster fique certinho dentro do Brasil
batsbrazil_forma <- mask(bats_brasraster, brazil_map)

plot(batsbrazil_forma)
plot(brazil_map, add=T)

Agora com uma resolução mais grossa (2 graus long-lat)

ras_brasil2 <- ras_brasil

res(ras_brasil2) <- 2

Raster de riqueza de morcegos

bats_brasraster2 <- rasterize(bats, ras_brasil2, field="binomial",fun=function(x,...){length(unique(na.omit(x)))})

Graficar

batsbrazil_forma2 <- mask(bats_brasraster2, brazil_map)

plot(batsbrazil_forma2)
plot(brazil_map, add=T)

Ficou legal? Mudou?

Resolução mais grossa ainda

ras_brasil5 <- ras_brasil2

res(ras_brasil5) <- 5

Raster de riqueza de morcegos

bats_brasraster5 <- rasterize(bats, ras_brasil5, field="binomial",fun=function(x,...){length(unique(na.omit(x)))})

Graficar.

batsbrazil_forma5 <- mask(bats_brasraster5, brazil_map)

plot(batsbrazil_forma5)
plot(brazil_map, add=T)

Mudou? Quanto? Aonde tem mais diferenças?

par(mfrow=c(1,3))

plot(batsbrazil_forma)
plot(brazil_map, add=T)

plot(batsbrazil_forma2)
plot(brazil_map, add=T)

plot(batsbrazil_forma5)
plot(brazil_map, add=T)

dev.off()