Análise Caatinga RL & APP
Felipe Melo & Adriana Pallegrini
04/06/2020
Usando os dados antigos ainda para consolidar resultados
Juntando todas as bases. São milhoes de propriedades privadas (PP). Mas como salvei as tabelas e dados finais no RData, não precisa correr essa parte do script
# filelist18 <- list.files(pattern = "*18.*.txt")
# datalist18 <- lapply(filelist18, function(x)read.table(x, header=T))
# PL_caat_2018 <- do.call("rbind", datalist18)
# PL_caat_2018$yr<-as.factor(rep(2018, each= n))
# filelist08 <- list.files(pattern = "*08.*.txt")
# datalist08 <- lapply(filelist08, function(x)read.table(x, header=T))
# PL_caat_2008 <- do.call("rbind", datalist18)
# PL_caat_2008$yr<-as.factor(rep(2008, each= n))
#
# data<-rbind(PL_caat_2008, PL_caat_2018)Aqui começa a brincadeira… Aqui estou criando várias variáveis que vamso usar para fazer os resumos por estado ou tamanho da PP. “Data” vai ser nossa tabela mãe da qual podemos extrair todas as informações que desejarmos. Importante checar se todas as descriçoes das variáveis fazem sentido bem como os summaries. Por exemplo, deve have ainda errinhos nos cálculos que acusam % maiores que 100% para “perc_caat_pp”, por exemplo
data %>%
mutate(AreaHa_caat2=AreaHa_caat*900/10000) %>%# transformando pixel em hectares da área que SERIA coberta com caatinga em cada PP segundo MapBiomas
mutate(vegHa=(X3+X4+X5+X11+X12+X13)*900/10000) %>% #cria a quant de vegetação em ha
mutate(dif_area=areaHa/AreaHa_caat2*100) %>% # Diferecença percentual entre área da propriedade e área de Caatinga. "Inf" quer dizer ZERO Caatinga segundo mapbimas
mutate(perc_veg_pp=vegHa/areaHa*100) %>% # cria a % de vegetação (incluidno outros biomas) por propriedade
mutate(perc_caat_pp=vegHa/AreaHa_caat2*100) %>% # cria a % de caatinga por propriedade. Aqui "NaN" é zero Caatinga
mutate_at(vars(perc_caat_pp), ~replace(., is.nan(.), 0)) %>% #tranforma NaN em zero
mutate(plSize2=as.factor(ifelse(grepl("S", plSize), "S", "L"))) %>% # transformar em grande e pequenos
mutate(vsd_cat=cut(perc_caat_pp, breaks=c(-Inf, 20, 40, 60, 80, Inf),
labels=c("less_20","20-40","40-60","60-80", "more_80"))) %>%
mutate(vsd_ha=areaHa*0.2) %>% # área que deveriam ter cada propriedade, 20%
mutate(net_vsd= vegHa - vsd_ha) %>% # Área de Caatinga que deveruia
mutate(dif_perc_veg=perc_veg_pp-perc_caat_pp) -> data # diferença entre vegetação já em porcentagem
summary(data) # podemos criar mais colunas com os cálculos que desejarmos## ID Count cd_mun
## Min. : 7000007 Min. : 1.0 Min. :2200053
## 1st Qu.: 26129932 1st Qu.: 11.0 1st Qu.:2312205
## Median :600628841 Median : 34.0 Median :2601805
## Mean :403887355 Mean : 261.3 Mean :2588957
## 3rd Qu.:601292854 3rd Qu.: 122.0 3rd Qu.:2901205
## Max. :605076304 Max. :469234.0 Max. :3171030
##
## nm_mun cd_micro nm_micro
## OURICURI : 41076 Min. :22001 ALTO MÉDIO CANINDÉ: 260236
## ARARIPINA : 37564 1st Qu.:23028 PAJEÚ : 197276
## MONTES CLAROS : 31456 Median :26001 ARARIPINA : 182184
## BRASÍLIA DE MINAS: 25436 Mean :25791 GARANHUNS : 126692
## PETROLINA : 25416 3rd Qu.:29002 GUANAMBI : 111208
## SANTA CRUZ : 25284 Max. :31007 MONTES CLAROS : 109412
## (Other) :5632292 (Other) :4831516
## cd_meso nm_meso cd_uf
## Min. :2201 SERTÃO PERNAMBUCANO : 489460 Min. :22.00
## 1st Qu.:2306 SUDESTE PIAUIENSE : 417088 1st Qu.:23.00
## Median :2601 AGRESTE PERNAMBUCANO: 409664 Median :26.00
## Mean :2581 CENTRO SUL BAIANO : 400764 Mean :25.78
## 3rd Qu.:2901 NORTE DE MINAS : 329800 3rd Qu.:29.00
## Max. :3102 SERTÃO PARAIBANO : 300480 Max. :31.00
## (Other) :3471268
## nm_uf nm_region biome country
## BA :1165192 NORDESTE:5488724 CAATINGA :4855684 BRASIL:5818524
## PE :1023272 SUDESTE : 329800 MATA ATLÂNTICA: 412620
## CE : 840176 CERRADO : 550220
## PI : 723480
## PB : 696584
## RN : 525744
## (Other): 844076
## landTenure pc_pl plSize areaHa
## AS : 8628 PL :5809624 L: 67920 Min. : 0.25
## CAR:1969180 PL1: 8900 M: 221496 1st Qu.: 2.75
## PR : 65572 S:5529108 Median : 8.50
## QL : 272 Mean : 65.32
## SI :3774872 3rd Qu.: 30.50
## Max. :117308.50
##
## areaMF X3 X4 X5
## Min. : 0.0028 Min. : 0.00 Min. : 0.0 Min. : 0.000
## 1st Qu.: 0.0571 1st Qu.: 0.00 1st Qu.: 0.0 1st Qu.: 0.000
## Median : 0.1750 Median : 0.00 Median : 0.0 Median : 0.000
## Mean : 1.2298 Mean : 18.31 Mean : 178.1 Mean : 0.069
## 3rd Qu.: 0.5962 3rd Qu.: 0.00 3rd Qu.: 0.0 3rd Qu.: 0.000
## Max. :1675.8357 Max. :197775.00 Max. :718209.0 Max. :10779.000
##
## X11 X12 X13 AreaHa_caat
## Min. :0 Min. : 0.0 Min. :0.00e+00 Min. : 0.0
## 1st Qu.:0 1st Qu.: 0.0 1st Qu.:0.00e+00 1st Qu.: 0.0
## Median :0 Median : 0.0 Median :0.00e+00 Median : 0.0
## Mean :0 Mean : 12.1 Mean :1.35e-03 Mean : 353.4
## 3rd Qu.:0 3rd Qu.: 0.0 3rd Qu.:0.00e+00 3rd Qu.: 0.0
## Max. :0 Max. :291032.0 Max. :3.01e+02 Max. :879492.0
##
## yr AreaHa_caat2 vegHa dif_area
## 2008:2909262 Min. : 0.0 Min. : 0.00 Min. :0.002619
## 2018:2909262 1st Qu.: 0.0 1st Qu.: 0.00 1st Qu.: Inf
## Median : 0.0 Median : 0.00 Median : Inf
## Mean : 31.8 Mean : 18.77 Mean : Inf
## 3rd Qu.: 0.0 3rd Qu.: 0.00 3rd Qu.: Inf
## Max. :79154.3 Max. :65605.23 Max. : Inf
##
## perc_veg_pp perc_caat_pp plSize2 dif_perc_veg
## Min. : 0.0 Min. : 0.000 L: 289416 Min. : -151.2
## 1st Qu.: 0.0 1st Qu.: 0.000 S:5529108 1st Qu.: 0.0
## Median : 0.0 Median : 0.000 Median : 0.0
## Mean : 178.1 Mean : 6.762 Mean : 171.4
## 3rd Qu.: 0.0 3rd Qu.: 0.000 3rd Qu.: 0.0
## Max. :2227236.0 Max. :156.977 Max. :2227159.2
##
## vsd_cat vsd_ha_GP vsd_ha vsd_cat_all
## less_20:5175076 Min. : 0.05 Min. : 0.05 Min. :-23461.70
## 20-40 : 177690 1st Qu.: 0.55 1st Qu.: 0.55 1st Qu.: -4.60
## 40-60 : 170816 Median : 1.70 Median : 1.70 Median : -1.05
## 60-80 : 163966 Mean : 13.06 Mean : 13.06 Mean : 18.74
## more_80: 130976 3rd Qu.: 6.10 3rd Qu.: 6.10 3rd Qu.: -0.20
## Max. :23461.70 Max. :23461.70 Max. : 79152.43
##
## net_vsd
## Min. :-23461.70
## 1st Qu.: -4.80
## Median : -1.20
## Mean : 5.71
## 3rd Qu.: -0.30
## Max. : 65170.74
## Agora vamos para algumas perguntas específicas:
Qual a % média de Caatinga por PP, segundo plSize e UF?
É curioso que a média percentual de Caatinga é constantemente menos em PP, o que as coloca em situação vulnerável em termos de serviços ecossistêmicos. Agora, a média é ainda abaixo do esperado em TODOS os estados paras as GP que são obrigadas ter 20% no mínimo.
## Perguntas específicas ----
# Qual a % média de Caatinga por PP, segundo plSize e UF?----
data %>%
filter(yr == "2018") %>% #Nesse caso devemos ficar com o ano de 2018
group_by(nm_uf, plSize2) %>%
summarise_at(vars(perc_caat_pp), funs(mean, se=sd(.)/sqrt(n()))) ->a # se correr até aqui e antes do %>%, dá pra ver a tabela## Warning: funs() is soft deprecated as of dplyr 0.8.0
## Please use a list of either functions or lambdas:
##
## # Simple named list:
## list(mean = mean, median = median)
##
## # Auto named with `tibble::lst()`:
## tibble::lst(mean, median)
##
## # Using lambdas
## list(~ mean(., trim = .2), ~ median(., na.rm = TRUE))
## This warning is displayed once per session.kable(a,align = "lccrr")| nm_uf | plSize2 | mean | se |
|---|---|---|---|
| AL | L | 7.199365 | 0.4646424 |
| AL | S | 5.257887 | 0.0422565 |
| PE | L | 10.192365 | 0.2193408 |
| PE | S | 9.223568 | 0.0295474 |
| BA | L | 7.284168 | 0.0916448 |
| BA | S | 3.451860 | 0.0191474 |
| MG | L | 12.657675 | 0.2886314 |
| MG | S | 6.403194 | 0.0403988 |
| SE | L | 7.391737 | 0.4019894 |
| SE | S | 1.397895 | 0.0214230 |
| CE | L | 11.112495 | 0.1696429 |
| CE | S | 6.057767 | 0.0318083 |
| RN | L | 18.992051 | 0.3145556 |
| RN | S | 4.498422 | 0.0310856 |
| PB | L | 13.037007 | 0.2817429 |
| PB | S | 6.757899 | 0.0341430 |
| PI | L | 12.803769 | 0.2523019 |
| PI | S | 11.714991 | 0.0453854 |
a %>%
ggplot(aes(nm_uf, mean, group=plSize2, fill=plSize2, color=plSize2))+
geom_col(position = "dodge")+
geom_errorbar(aes(ymin=mean-se, ymax=mean+se),position=position_dodge(width=0.9), colour="black", width=.2)+
labs(x="States with Caatinga", y = "Mean % of Caatinga per private property")+
theme(axis.text.x=element_blank(), axis.title.x=element_blank())-> graf1
graf1
Qual a quantidade de Caatinga por PP, segundo plSize e UF?
data %>%
filter(yr == "2018") %>% #Nesse caso devemos ficar com o ano de 2018
group_by(nm_uf, plSize2) %>%
summarise_at(vars(vegHa), funs(sum)) %>% # se correr até aqui e antes do %>%, dá pra ver a tabela
ggplot(aes(nm_uf, vegHa, group=plSize2, fill=plSize2, color=plSize2))+
geom_col(position = "dodge")+
labs(x="States with Caatinga", y = "Amount of Caatinga vegetation (ha)") -> graf2
graf2
#Juntando os dois gráfico para ter uma idea
plot_grid(graf1,graf2, nrow = 2, align = TRUE) # Aqui podemos ter uma idea de como é a distribuição das florestas por estados
Qual a distribuição do número de PP segundo % de Caatinga por plSize e UF?
data %>%
filter(yr == "2018") %>%
group_by(nm_uf, plSize2,vsd_cat) %>%
summarise(count=n()) %>%
ggplot(aes(nm_uf, log(count), group=vsd_cat, fill=vsd_cat, color=vsd_cat))+
geom_col(position = "dodge")+
facet_wrap(.~plSize2)
Aqui em tabela
data %>%
filter(yr == "2018") %>%
group_by(nm_uf,plSize2,vsd_cat) %>%
summarise(count=n()) %>%
spread(vsd_cat,count) %>%
kable()| nm_uf | plSize2 | less_20 | 20-40 | 40-60 | 60-80 | more_80 |
|---|---|---|---|---|---|---|
| AL | L | 1352 | 100 | 52 | 40 | 20 |
| AL | S | 130147 | 5566 | 3481 | 2232 | 1718 |
| PE | L | 11260 | 313 | 379 | 572 | 706 |
| PE | S | 416705 | 26832 | 25689 | 20965 | 8215 |
| BA | L | 53461 | 1215 | 1244 | 1211 | 2589 |
| BA | S | 491802 | 11878 | 8176 | 5681 | 5339 |
| MG | L | 5453 | 582 | 456 | 316 | 231 |
| MG | S | 138835 | 10244 | 5036 | 2289 | 1458 |
| SE | L | 1390 | 116 | 72 | 36 | 6 |
| SE | S | 108075 | 1636 | 692 | 303 | 104 |
| CE | L | 23288 | 286 | 490 | 952 | 2174 |
| CE | S | 357505 | 4654 | 7825 | 12611 | 10303 |
| RN | L | 6854 | 462 | 684 | 916 | 724 |
| RN | S | 233275 | 6400 | 6083 | 4754 | 2720 |
| PB | L | 7375 | 402 | 370 | 582 | 505 |
| PB | S | 302008 | 9661 | 9648 | 9523 | 8218 |
| PI | L | 13198 | 34 | 140 | 389 | 1711 |
| PI | S | 285555 | 8464 | 14891 | 18611 | 18747 |
Quanta Caatinga falta por UF separado por PlSize2
Aqui estamos usando o superavit de pequenas para fazer
Nem todos os estados tem deficit ou todas as categorias de PP (L and S) têm déficit
data %>%
filter(yr == "2018") %>%
group_by(nm_uf,vsd_cat, plSize2) %>%
summarise_at(vars(net_vsd), funs(sum)) %>%
spread(vsd_cat, net_vsd) %>%
as.data.frame() %>%
mutate(VSD_GP = rowSums(.[4:7])+ rowSums(.[3])) %>%
kable()| nm_uf | plSize2 | less_20 | 20-40 | 40-60 | 60-80 | more_80 | VSD_GP |
|---|---|---|---|---|---|---|---|
| AL | L | -160894.6 | 5605.88 | 11290.22 | 11239.81 | 6207.90 | -126550.81 |
| AL | S | -256455.1 | 80521.22 | 62557.58 | 19531.24 | 13905.63 | -79939.41 |
| PE | L | -1134038.2 | 39801.26 | 119792.57 | 331346.22 | 448992.35 | -194105.77 |
| PE | S | -1287636.1 | 320680.09 | 852204.96 | 1668005.84 | 961300.56 | 2514555.31 |
| BA | L | -10276505.3 | 464076.71 | 1192263.71 | 1661385.84 | 4465274.94 | -2493504.13 |
| BA | S | -3140755.5 | 1577262.61 | 2809306.07 | 2964252.99 | 2517382.91 | 6727449.12 |
| MG | L | -920176.7 | 47007.34 | 158717.01 | 167571.15 | 97328.57 | -449552.59 |
| MG | S | -551318.7 | 76604.35 | 155968.12 | 145568.82 | 86062.74 | -87114.70 |
| SE | L | -131916.6 | 16091.74 | 15779.56 | 12552.20 | 1403.06 | -86090.02 |
| SE | S | -270913.3 | 61126.63 | 43975.56 | 31630.46 | 7385.42 | -126795.21 |
| CE | L | -2813236.5 | 61568.71 | 247526.07 | 808823.73 | 2185169.62 | 489851.59 |
| CE | S | -1682600.0 | 140568.53 | 663070.88 | 2238024.88 | 2892194.06 | 4251258.31 |
| RN | L | -651687.2 | 48816.45 | 229262.93 | 499438.47 | 462192.70 | 588023.31 |
| RN | S | -595343.6 | 192430.40 | 430851.80 | 562406.95 | 381898.14 | 972243.66 |
| PB | L | -754692.9 | 35128.72 | 90605.57 | 249429.68 | 263918.00 | -115610.92 |
| PB | S | -823782.7 | 224643.29 | 514071.59 | 884938.53 | 715329.82 | 1515200.52 |
| PI | L | -2806458.2 | 4333.27 | 55620.73 | 428356.63 | 2237317.23 | -80830.37 |
| PI | S | -1718168.0 | 67777.42 | 360392.37 | 1245988.96 | 3429117.38 | 3385108.10 |
Num plot se vê melhor
data %>%
filter(yr == "2018") %>%
group_by(nm_uf,vsd_cat, plSize2) %>%
summarise_at(vars(net_vsd), funs(sum)) %>%
spread(vsd_cat, net_vsd) %>%
as.data.frame() %>%
mutate(VSD_GP = rowSums(.[4:7])+ rowSums(.[3])) %>% # Nem todos os estados tem deficit ou todas as categorias de PP (L and S) têm déficit
ggplot(aes(nm_uf,VSD_GP, group=plSize2, fill=plSize2))+
geom_col(position = "dodge")+coord_flip()
Quanta Caatinga falta na categoria <20% por UF separada por PlSize2
data %>%
filter(yr == "2018") %>%
group_by(nm_uf,vsd_cat, plSize2) %>%
summarise_at(vars(net_vsd), funs(sum)) %>%
spread(vsd_cat, net_vsd) %>%
as.data.frame() %>%
ggplot(aes(nm_uf,less_20, group=plSize2, fill=plSize2))+
geom_col(position = "dodge")+coord_flip()
Quais tamanhos contribuem mais com o déficit de Caatinga?
data %>%
filter(yr == "2018") %>%
group_by(vsd_cat, plSize2) %>%
summarise_at(vars(net_vsd), funs(sum)) %>%
spread(vsd_cat, net_vsd) %>%
mutate(VSD_GP = rowSums(.[3:6])+ rowSums(.[2])) %>%
kable()| plSize2 | less_20 | 20-40 | 40-60 | 60-80 | more_80 | VSD_GP |
|---|---|---|---|---|---|---|
| L | -19649606 | 722430.1 | 2120858 | 4170144 | 10167804 | -2468370 |
| S | -10326973 | 2741614.5 | 5892399 | 9760349 | 11004577 | 19071966 |