Resultados Modelos

library(tidyverse)
library(lme4)
library(parameters)

rsultado1 <- readRDS("~/RStudio/ENEM/Resultados/Resultado_idh_prop_1_rep_1_modelo_bs_sex.rds")
rsultado2 <- readRDS("~/RStudio/ENEM/Resultados/Resultado_idh_Tipo_Ruralxurbano_prop_1_rep_1_modelo_bs_sex.rds")
rsultado3 <- readRDS("~/RStudio/ENEM/Resultados/Resultado_idh_Tipo_Escola_prop_1_rep_1_modelo_bs_sex.rds")
rsultado4 <- readRDS("~/RStudio/ENEM/Resultados/Resultado_idh_Tipo_Sexo_prop_1_rep_1_modelo_bs_sex.rds")
rsultado5 <- readRDS("~/RStudio/ENEM/Resultados/Resultado_idh_TP_2009_prop_1_rep_1_modelo_bs_sex.rds")

Modelo Basico

resultado<-rsultado1
resultado$formula[1]#Nota=Poluente+Renda+educaçãoparental+ANO+Sexo+Loc_escola+(1+poluente|Estado)

[1] "samp[[y]] ~ samp[[x]] + Q15 + Q10 + IDHM + TP_DEPENDENCIA_ADM_ESC + bs(as.numeric(NU_ANO)) + TP_SEXO + TP_LOCALIZACAO_ESC + (1 + samp[[x]] | SG_UF_RESIDENCIA)"

summary(as.factor(resultado$amostragem))

15443772 
     162 

R <-resultado %>% mutate(lowerCI=slope-1.96*SE_Slope_rand,upperCI=slope+1.96*SE_Slope_rand) %>% 
  mutate(sinal=ifelse(lowerCI>0,"Positivo",ifelse(upperCI<0,"Negativo","null")))
#print(table(R$poluente,R$sinal))
table(R$modelo,R$sinal)

                            
                             Negativo null Positivo
  no2_ppb.NU_NOTA_OBJETIVA         10    7       10
  no2_ppb.NU_NOTA_REDACAO          13    6        8
  o3_ppb.NU_NOTA_OBJETIVA          22    4        1
  o3_ppb.NU_NOTA_REDACAO           20    5        2
  pm25_ugm3.NU_NOTA_OBJETIVA       11    8        8
  pm25_ugm3.NU_NOTA_REDACAO        10    9        8

R %>%   ggplot(aes(x=slope,y=fator,xmin=lowerCI,xmax=upperCI,col=sinal))+
  #geom_point(alpha=0.5)+
  geom_pointrange(alpha=0.5)+facet_grid(nota~poluente,scales = "free")+
  geom_vline(xintercept = 0, color = "red", size=0.1,alpha=5)+
  #geom_errorbar(alpha=0.5)+
  labs(title = "Modelo mistos",
       subtitle = paste0("n=",R$amostragem,"
",R$formula),
       caption = "EPPG FGV")+xlab("Slope")+ylab("Estado")

resultado$sinal<-ifelse(resultado$lowerCI_efx>0,"Positivo",ifelse(resultado$upperCI_efx<0,"Negativo","null"))
R2<-resultado[6:19] %>% distinct() 
table(R2$modelo,R2$sinal)

                            
                             Negativo null
  no2_ppb.NU_NOTA_OBJETIVA          0    1
  no2_ppb.NU_NOTA_REDACAO           0    1
  o3_ppb.NU_NOTA_OBJETIVA           1    0
  o3_ppb.NU_NOTA_REDACAO            1    0
  pm25_ugm3.NU_NOTA_OBJETIVA        0    1
  pm25_ugm3.NU_NOTA_REDACAO         0    1

R2 %>%   ggplot(aes(x=efeito_fixo,y=teste,xmin=lowerCI_efx,xmax=upperCI_efx,col=sinal))+geom_point(alpha=0.5)+
  geom_pointrange(alpha=0.5)+facet_grid(nota~poluente,scales = "free")+
  geom_vline(xintercept = 0, color = "red", size=0.1,alpha=0.5)+
  #geom_errorbar(alpha=0.5)+
  labs(title = "Modelo mistos",
       subtitle = paste0("n=",R$amostragem,"
",R$formula),
       caption = "EPPG FGV")+xlab("Slope")+ylab("Estado")

Modelo com Subset Rural Urbano

resultado<-rsultado2
resultado$formula[1]#Nota=Poluente+Renda+educaçãoparental+ANO+Sexo+Loc_escola+(1+poluente|Estado)

[1] "samp[[y]] ~ samp[[x]] + Q15 + Q10 + IDHM + TP_DEPENDENCIA_ADM_ESC + bs(as.numeric(NU_ANO)) + TP_SEXO + (1 + samp[[x]] | SG_UF_RESIDENCIA)"

summary(as.factor(resultado$amostragem))

 178002 8836785 
    156     162 

R <-resultado %>% mutate(lowerCI=slope-1.96*SE_Slope_rand,upperCI=slope+1.96*SE_Slope_rand) %>% 
  mutate(sinal=ifelse(lowerCI>0,"Positivo",ifelse(upperCI<0,"Negativo","null")))
#print(table(R$poluente,R$sinal))
table(R$modelo,R$sinal)

                            
                             Negativo null Positivo
  no2_ppb.NU_NOTA_OBJETIVA         21   14       18
  no2_ppb.NU_NOTA_REDACAO          19   17       17
  o3_ppb.NU_NOTA_OBJETIVA          40   11        2
  o3_ppb.NU_NOTA_REDACAO           37   12        4
  pm25_ugm3.NU_NOTA_OBJETIVA       23   14       16
  pm25_ugm3.NU_NOTA_REDACAO        22   18       13

table(R$subset,R$sinal,R$poluente)

, ,  = no2_ppb

   
    Negativo null Positivo
  1       25   11       18
  2       15   20       17

, ,  = o3_ppb

   
    Negativo null Positivo
  1       42    9        3
  2       35   14        3

, ,  = pm25_ugm3

   
    Negativo null Positivo
  1       21   16       17
  2       24   16       12

R %>%   ggplot(aes(x=slope,y=fator,xmin=lowerCI,xmax=upperCI,col=sinal))+
  #geom_point(alpha=0.5)+
  geom_pointrange(alpha=0.5)+facet_grid(nota*subset~poluente,scales = "free")+
  geom_vline(xintercept = 0, color = "red", size=0.1,alpha=5)+
  #geom_errorbar(alpha=0.5)+
  labs(title = "Modelo mistos",
       subtitle = paste0("n=",R$amostragem,"
",R$formula),
       caption = "EPPG FGV")+xlab("Slope")+ylab("Estado")

resultado$sinal<-ifelse(resultado$lowerCI_efx>0,"Positivo",ifelse(resultado$upperCI_efx<0,"Negativo","null"))
R2<-resultado[6:19] %>% distinct() 
table(R2$modelo,R2$sinal)

                            
                             Negativo null
  no2_ppb.NU_NOTA_OBJETIVA          0    2
  no2_ppb.NU_NOTA_REDACAO           0    2
  o3_ppb.NU_NOTA_OBJETIVA           2    0
  o3_ppb.NU_NOTA_REDACAO            2    0
  pm25_ugm3.NU_NOTA_OBJETIVA        0    2
  pm25_ugm3.NU_NOTA_REDACAO         0    2

R2 %>%   ggplot(aes(x=efeito_fixo,y=teste,xmin=lowerCI_efx,xmax=upperCI_efx,col=sinal))+geom_point(alpha=0.5)+
  geom_pointrange(alpha=0.5)+facet_grid(nota*subset~poluente,scales = "free")+
  geom_vline(xintercept = 0, color = "red", size=0.1,alpha=0.5)+
  #geom_errorbar(alpha=0.5)+
  labs(title = "Modelo mistos",
       subtitle = paste0("n=",R$amostragem,"
",R$formula),
       caption = "EPPG FGV")+xlab("Slope")+ylab("Estado")

Modelo com subset Tipo de Escola

resultado<-rsultado3
resultado$formula[1]#Nota=Poluente+Renda+educaçãoparental+ANO+Sexo+Loc_escola+(1+poluente|Estado)

[1] "samp[[y]] ~ samp[[x]] + Q15 + Q10 + IDHM + bs(as.numeric(NU_ANO)) + TP_SEXO + TP_LOCALIZACAO_ESC + (1 + samp[[x]] | SG_UF_RESIDENCIA)"

summary(as.factor(resultado$amostragem))

  239813   302601  2825635 12075723 
     150      138      162      162 

R <-resultado %>% mutate(lowerCI=slope-1.96*SE_Slope_rand,upperCI=slope+1.96*SE_Slope_rand) %>% 
  mutate(sinal=ifelse(lowerCI>0,"Positivo",ifelse(upperCI<0,"Negativo","null")))
#print(table(R$poluente,R$sinal))
table(R$modelo,R$sinal)

                            
                             Negativo null Positivo
  no2_ppb.NU_NOTA_OBJETIVA         35   30       37
  no2_ppb.NU_NOTA_REDACAO          35   34       33
  o3_ppb.NU_NOTA_OBJETIVA          76   22        4
  o3_ppb.NU_NOTA_REDACAO           67   27        8
  pm25_ugm3.NU_NOTA_OBJETIVA       39   29       34
  pm25_ugm3.NU_NOTA_REDACAO        42   31       29

table(R$subset,R$sinal,R$poluente)

, ,  = no2_ppb

   
    Negativo null Positivo
  1       12   20       18
  2       24   13       17
  3        9   17       20
  4       25   14       15

, ,  = o3_ppb

   
    Negativo null Positivo
  1       34   13        3
  2       42   10        2
  3       28   14        4
  4       39   12        3

, ,  = pm25_ugm3

   
    Negativo null Positivo
  1       16   18       16
  2       23   14       17
  3       16   14       16
  4       26   14       14

R %>%   ggplot(aes(x=slope,y=fator,xmin=lowerCI,xmax=upperCI,col=sinal))+
  #geom_point(alpha=0.5)+
  geom_pointrange(alpha=0.5)+facet_grid(nota*subset~poluente,scales = "free")+
  geom_vline(xintercept = 0, color = "red", size=0.1,alpha=5)+
  #geom_errorbar(alpha=0.5)+
  labs(title = "Modelo mistos",
       subtitle = paste0("n=",R$amostragem,"
",R$formula),
       caption = "EPPG FGV")+xlab("Slope")+ylab("Estado")

resultado$sinal<-ifelse(resultado$lowerCI_efx>0,"Positivo",ifelse(resultado$upperCI_efx<0,"Negativo","null"))
R2<-resultado[6:19] %>% distinct() 
table(R2$modelo,R2$sinal)

                            
                             Negativo null
  no2_ppb.NU_NOTA_OBJETIVA          0    4
  no2_ppb.NU_NOTA_REDACAO           0    4
  o3_ppb.NU_NOTA_OBJETIVA           4    0
  o3_ppb.NU_NOTA_REDACAO            4    0
  pm25_ugm3.NU_NOTA_OBJETIVA        0    4
  pm25_ugm3.NU_NOTA_REDACAO         0    4

R2 %>%   ggplot(aes(x=efeito_fixo,y=teste,xmin=lowerCI_efx,xmax=upperCI_efx,col=sinal))+geom_point(alpha=0.5)+
  geom_pointrange(alpha=0.5)+facet_grid(nota*subset~poluente,scales = "free")+
  geom_vline(xintercept = 0, color = "red", size=0.1,alpha=0.5)+
  #geom_errorbar(alpha=0.5)+
  labs(title = "Modelo mistos",
       subtitle = paste0("n=",R$amostragem,"
",R$formula),
       caption = "EPPG FGV")+xlab("Slope")+ylab("Estado")

Modelo com subset Sexo

resultado<-rsultado4
resultado$formula[1]#Nota=Poluente+Renda+educaçãoparental+ANO+Sexo+Loc_escola+(1+poluente|Estado)

[1] "samp[[y]] ~ samp[[x]] + Q15 + Q10 + IDHM + TP_DEPENDENCIA_ADM_ESC + bs(as.numeric(NU_ANO)) + +TP_LOCALIZACAO_ESC + (1 + samp[[x]] | SG_UF_RESIDENCIA)"

summary(as.factor(resultado$amostragem))

3983927 5030860 
    162     162 

R <-resultado %>% mutate(lowerCI=slope-1.96*SE_Slope_rand,upperCI=slope+1.96*SE_Slope_rand) %>% 
  mutate(sinal=ifelse(lowerCI>0,"Positivo",ifelse(upperCI<0,"Negativo","null")))
#print(table(R$poluente,R$sinal))
table(R$modelo,R$sinal)

                            
                             Negativo null Positivo
  no2_ppb.NU_NOTA_OBJETIVA         21   13       20
  no2_ppb.NU_NOTA_REDACAO          25   13       16
  o3_ppb.NU_NOTA_OBJETIVA          42   10        2
  o3_ppb.NU_NOTA_REDACAO           38   12        4
  pm25_ugm3.NU_NOTA_OBJETIVA       21   19       14
  pm25_ugm3.NU_NOTA_REDACAO        20   20       14

table(R$subset,R$sinal,R$poluente)

, ,  = no2_ppb

   
    Negativo null Positivo
  F       22   14       18
  M       24   12       18

, ,  = o3_ppb

   
    Negativo null Positivo
  F       39   12        3
  M       41   10        3

, ,  = pm25_ugm3

   
    Negativo null Positivo
  F       22   19       13
  M       19   20       15

R %>%   ggplot(aes(x=slope,y=fator,xmin=lowerCI,xmax=upperCI,col=sinal))+
  #geom_point(alpha=0.5)+
  geom_pointrange(alpha=0.5)+facet_grid(nota*subset~poluente,scales = "free")+
  geom_vline(xintercept = 0, color = "red", size=0.1,alpha=5)+
  #geom_errorbar(alpha=0.5)+
  labs(title = "Modelo mistos",
       subtitle = paste0("n=",R$amostragem,"
",R$formula),
       caption = "EPPG FGV")+xlab("Slope")+ylab("Estado")

resultado$sinal<-ifelse(resultado$lowerCI_efx>0,"Positivo",ifelse(resultado$upperCI_efx<0,"Negativo","null"))
R2<-resultado[6:19] %>% distinct() 
table(R2$modelo,R2$sinal)

                            
                             Negativo null
  no2_ppb.NU_NOTA_OBJETIVA          0    2
  no2_ppb.NU_NOTA_REDACAO           0    2
  o3_ppb.NU_NOTA_OBJETIVA           2    0
  o3_ppb.NU_NOTA_REDACAO            2    0
  pm25_ugm3.NU_NOTA_OBJETIVA        0    2
  pm25_ugm3.NU_NOTA_REDACAO         0    2

R2 %>%   ggplot(aes(x=efeito_fixo,y=teste,xmin=lowerCI_efx,xmax=upperCI_efx,col=sinal))+geom_point(alpha=0.5)+
  geom_pointrange(alpha=0.5)+facet_grid(nota*subset~poluente,scales = "free")+
  geom_vline(xintercept = 0, color = "red", size=0.1,alpha=0.5)+
  #geom_errorbar(alpha=0.5)+
  labs(title = "Modelo mistos",
       subtitle = paste0("n=",R$amostragem,"
",R$formula),
       caption = "EPPG FGV")+xlab("Slope")+ylab("Estado")

Modelo com subset por Renda

resultado<-rsultado5
resultado$formula[1]#Nota=Poluente+Renda+educaçãoparental+ANO+Sexo+Loc_escola+(1+poluente|Estado)

[1] "samp[[y]] ~ samp[[x]] + Q15 + Q10 + IDHM + TP_DEPENDENCIA_ADM_ESC + bs(as.numeric(NU_ANO)) + TP_SEXO + TP_LOCALIZACAO_ESC + (1 + samp[[x]] | SG_UF_RESIDENCIA)"

summary(as.factor(resultado$amostragem))

2581321 6433466 
    126     162 

R <-resultado %>% mutate(lowerCI=slope-1.96*SE_Slope_rand,upperCI=slope+1.96*SE_Slope_rand) %>% 
  mutate(sinal=ifelse(lowerCI>0,"Positivo",ifelse(upperCI<0,"Negativo","null")))
#print(table(R$poluente,R$sinal))
table(R$modelo,R$sinal)

                            
                             Negativo null Positivo
  no2_ppb.NU_NOTA_OBJETIVA         14   12       22
  no2_ppb.NU_NOTA_REDACAO          16   18       14
  o3_ppb.NU_NOTA_OBJETIVA          25    7       16
  o3_ppb.NU_NOTA_REDACAO           25    7       16
  pm25_ugm3.NU_NOTA_OBJETIVA       19    8       21
  pm25_ugm3.NU_NOTA_REDACAO        23    9       16

table(R$subset,R$sinal,R$poluente)

, ,  = no2_ppb

          
           Negativo null Positivo
  pós-2008       14   16       24
  pré=2008       16   14       12

, ,  = o3_ppb

          
           Negativo null Positivo
  pós-2008       33    7       14
  pré=2008       17    7       18

, ,  = pm25_ugm3

          
           Negativo null Positivo
  pós-2008       18   10       26
  pré=2008       24    7       11

R %>%   ggplot(aes(x=slope,y=fator,xmin=lowerCI,xmax=upperCI,col=sinal))+
  #geom_point(alpha=0.5)+
  geom_pointrange(alpha=0.5)+facet_grid(nota*subset~poluente,scales = "free")+
  geom_vline(xintercept = 0, color = "red", size=0.1,alpha=5)+
  #geom_errorbar(alpha=0.5)+
  labs(title = "Modelo mistos",
       subtitle = paste0("n=",R$amostragem,"
",R$formula),
       caption = "EPPG FGV")+xlab("Slope")+ylab("Estado")

resultado$sinal<-ifelse(resultado$lowerCI_efx>0,"Positivo",ifelse(resultado$upperCI_efx<0,"Negativo","null"))
R2<-resultado[6:19] %>% distinct() 
table(R2$modelo,R2$sinal)

                            
                             Negativo null
  no2_ppb.NU_NOTA_OBJETIVA          0    2
  no2_ppb.NU_NOTA_REDACAO           0    2
  o3_ppb.NU_NOTA_OBJETIVA           1    1
  o3_ppb.NU_NOTA_REDACAO            1    1
  pm25_ugm3.NU_NOTA_OBJETIVA        0    2
  pm25_ugm3.NU_NOTA_REDACAO         0    2

R2 %>%   ggplot(aes(x=efeito_fixo,y=teste,xmin=lowerCI_efx,xmax=upperCI_efx,col=sinal))+geom_point(alpha=0.5)+
  geom_pointrange(alpha=0.5)+facet_grid(nota*subset~poluente,scales = "free")+
  geom_vline(xintercept = 0, color = "red", size=0.1,alpha=0.5)+
  #geom_errorbar(alpha=0.5)+
  labs(title = "Modelo mistos",
       subtitle = paste0("n=",R$amostragem,"
",R$formula),
       caption = "EPPG FGV")+xlab("Slope")+ylab("Estado")

Modelo com subset por Grupo de Anos

Observações:

• Padronizar o tamanho da amostra para análises: aumentar(tentar chegar a 100%)

• Esquecer repetições: colocar só um modelo nos graficos

• Atrito: quantos retirou por causa dos NA’s

• incluir variavel de Renda ou IDH no modelo