Dicionário de variáveis: https://basedosdados.org/dataset/86bac6cc-575f-4289-a857-13f3f52c9a1d?table=f8b6030c-3fb1-4f64-81ef-3c5c3c888b6c
Base$variables <- transform(Base$variables, Idade=C008)
Base$variables <- transform(Base$variables, Sexo=C006)
Base$variables <- transform(Base$variables, HorasSemanaisTrabalho=E017)
Base$variables <- transform(Base$variables, Raça=C009)
Base$variables <- transform(Base$variables, Rendimento_bruto=E01602)
Base$variables <- transform(Base$variables, Doença=J007)
Base$variables <- transform(Base$variables, Motivo_Saude=J00402)
Base$variables <- transform(Base$variables, Interrupção_Atividade=J002)
Base$variables <- transform(Base$variables, Nível_Instrução=VDD004A)
Base$variables <- transform(Base$variables, TemPlano=I00102)limpar_coluna_binaria <- function(base, nome_coluna) {
vars <- base$variables
# Converte para texto (remove restrição de fator)
vars[[nome_coluna]] <- as.character(vars[[nome_coluna]])
# Substitui "Sim" -> 1 e "Não" -> 0
vars[[nome_coluna]][vars[[nome_coluna]] == "Sim"] <- 1
vars[[nome_coluna]][vars[[nome_coluna]] == "Não"] <- 0
# Converte para numérico
vars[[nome_coluna]] <- as.numeric(vars[[nome_coluna]])
# Filtra apenas 0 e 1
vars <- vars[vars[[nome_coluna]] %in% c(0, 1), ]
# Atualiza o objeto survey.design
base$variables <- vars
base <- subset(base, vars[[nome_coluna]] %in% c(0, 1))
return(base)
}
Base = limpar_coluna_binaria(Base, "Doença")J002: Interrupção na realização de quaisquer atividades habituais por motivo da própria saúde
J00402: Motivo de saúde que impediu a realização das atividades habituais
## Total amostral:
#Variáveis Numéricas:
total_renda = svytotal(x=~Rendimento_bruto, design=Base, na.rm=TRUE)
total_renda total SE
Rendimento_bruto 2.1323e+11 3604489769
2.5 % 97.5 %
Rendimento_bruto 206163204561 220292544822
#Variáveis categóricas - Utilizando mais de uma coluna
svytotal(x=~Raça+Sexo, design=Base, na.rm=TRUE) total SE
RaçaBranca 91037722 686604.4
RaçaPreta 21786515 338021.3
RaçaAmarela 1674691 144663.7
RaçaParda 94086142 600647.8
RaçaIndígena 990265 74430.0
RaçaIgnorado 14272 5707.5
SexoHomem 100158109 240473.3
SexoMulher 109431498 240473.3
#Variáveis categóricas - Estimação por cruzamento
svytotal(x=~interaction(Sexo,Raça), design=Base, na.rm=TRUE) total SE
interaction(Sexo, Raça)Homem.Branca 42682904.8 378164.8
interaction(Sexo, Raça)Mulher.Branca 48354816.8 404879.1
interaction(Sexo, Raça)Homem.Preta 10691164.0 186451.9
interaction(Sexo, Raça)Mulher.Preta 11095350.7 200050.6
interaction(Sexo, Raça)Homem.Amarela 754080.8 73347.6
interaction(Sexo, Raça)Mulher.Amarela 920610.6 81402.9
interaction(Sexo, Raça)Homem.Parda 45496265.5 337575.6
interaction(Sexo, Raça)Mulher.Parda 48589876.7 351852.0
interaction(Sexo, Raça)Homem.Indígena 523966.4 52076.8
interaction(Sexo, Raça)Mulher.Indígena 466298.3 33901.4
interaction(Sexo, Raça)Homem.Ignorado 9727.6 3936.5
interaction(Sexo, Raça)Mulher.Ignorado 4544.8 1978.7
mean SE
Rendimento_bruto 2224.8 35.465
mean SE
SexoHomem 0.47788 0.0011
SexoMulher 0.52212 0.0011
variance SE
Rendimento_bruto 14048208 1720046
##Razão
#svyratio()
razão <- svyratio(numerator=~(Motivo_Saude=="Problemas respiratórios (Resfriado / gripe /sinusite/ asma / bronquite / pneumonia)"),
denominator=~(Interrupção_Atividade=="Sim"), design=Base, na.rm=TRUE)
razãoRatio estimator: svyratio.survey.design2(numerator = ~(Motivo_Saude == "Problemas respiratórios (Resfriado / gripe /sinusite/ asma / bronquite / pneumonia)"),
denominator = ~(Interrupção_Atividade == "Sim"), design = Base,
na.rm = TRUE)
Ratios=
Interrupção_Atividade == "Sim"
Motivo_Saude == "Problemas respiratórios (Resfriado / gripe /sinusite/ asma / bronquite / pneumonia)" 0.2096767
SEs=
Interrupção_Atividade == "Sim"
Motivo_Saude == "Problemas respiratórios (Resfriado / gripe /sinusite/ asma / bronquite / pneumonia)" 0.005052547
Interrupção_Atividade == "Sim"
Motivo_Saude == "Problemas respiratórios (Resfriado / gripe /sinusite/ asma / bronquite / pneumonia)" 0.02409685
2.5 %
Motivo_Saude == "Problemas respiratórios (Resfriado / gripe /sinusite/ asma / bronquite / pneumonia)"/Interrupção_Atividade == "Sim" 0.1997739
97.5 %
Motivo_Saude == "Problemas respiratórios (Resfriado / gripe /sinusite/ asma / bronquite / pneumonia)"/Interrupção_Atividade == "Sim" 0.2195795
Para obtenção do intervalo de confiança pela função svyquantile, é necessário utilizar o parâmetro ci=TRUE e a partir dele, utilizar as funções SE e cv para estimar o erro padrão e coeficiente de variação, respectivamente.
mediana = svyquantile(x=~Rendimento_bruto, design=Base, quantiles=0.5, ci=TRUE, na.rm=TRUE) ##MEDIANA
mediana$Rendimento_bruto
quantile ci.2.5 ci.97.5 se
0.5 1350 1300 1400 25.50653
attr(,"hasci")
[1] TRUE
attr(,"class")
[1] "newsvyquantile"
Rendimento_bruto
25.50653
Rendimento_bruto
0.01889373
svyquantile(x=~Rendimento_bruto, design=Base, quantiles=c(.25,.5,.75), ci=FALSE, na.rm=TRUE) #Múltiplos quantis$Rendimento_bruto
0.25 0.5 0.75
[1,] 998 1350 2200
attr(,"hasci")
[1] FALSE
attr(,"class")
[1] "newsvyquantile"
# Estimação por dominio (filtro)
svytotal(x=~Rendimento_bruto, design=subset(Base, Sexo=="Mulher"), na.rm=TRUE) total SE
Rendimento_bruto 8.1577e+10 1440084928
svytotal(x=~Nível_Instrução, design=subset(Base, Sexo=="Homem" & Raça=="Parda" & Idade>30), na.rm=TRUE) total SE
Nível_InstruçãoSem instrução 2573476 63547
Nível_InstruçãoFundamental incompleto ou equivalente 9167518 129978
Nível_InstruçãoFundamental completo ou equivalente 1943169 64405
Nível_InstruçãoMédio incompleto ou equivalente 1262245 48120
Nível_InstruçãoMédio completo ou equivalente 5777822 107139
Nível_InstruçãoSuperior incompleto ou equivalente 522290 30839
Nível_InstruçãoSuperior completo 2011870 62089
svymean(x=~Nível_Instrução, design=subset(Base, Sexo=="Homem" & Raça=="Parda" & Idade>30), na.rm=TRUE) mean SE
Nível_InstruçãoSem instrução 0.110647 0.0026
Nível_InstruçãoFundamental incompleto ou equivalente 0.394160 0.0045
Nível_InstruçãoFundamental completo ou equivalente 0.083547 0.0027
Nível_InstruçãoMédio incompleto ou equivalente 0.054271 0.0020
Nível_InstruçãoMédio completo ou equivalente 0.248419 0.0040
Nível_InstruçãoSuperior incompleto ou equivalente 0.022456 0.0013
Nível_InstruçãoSuperior completo 0.086501 0.0026
##Coeficiente de variação
cv(svytotal(x=~Rendimento_bruto, design=subset(Base, Sexo=="Mulher"), na.rm=TRUE)) Rendimento_bruto
Rendimento_bruto 0.01765309
#Estimação para vários dominios ("Groupby/Sumarize")
svyby(formula=~Sexo, by=~Nível_Instrução, design=Base, FUN=svymean, na.rm=TRUE) Nível_Instrução
Sem instrução Sem instrução
Fundamental incompleto ou equivalente Fundamental incompleto ou equivalente
Fundamental completo ou equivalente Fundamental completo ou equivalente
Médio incompleto ou equivalente Médio incompleto ou equivalente
Médio completo ou equivalente Médio completo ou equivalente
Superior incompleto ou equivalente Superior incompleto ou equivalente
Superior completo Superior completo
SexoHomem SexoMulher se.SexoHomem
Sem instrução 0.4892280 0.5107720 0.004866419
Fundamental incompleto ou equivalente 0.5013440 0.4986560 0.002143222
Fundamental completo ou equivalente 0.4919528 0.5080472 0.005334081
Médio incompleto ou equivalente 0.5078484 0.4921516 0.005893286
Médio completo ou equivalente 0.4636716 0.5363284 0.003156167
Superior incompleto ou equivalente 0.4635797 0.5364203 0.007745619
Superior completo 0.4036566 0.5963434 0.003870122
se.SexoMulher
Sem instrução 0.004866419
Fundamental incompleto ou equivalente 0.002143222
Fundamental completo ou equivalente 0.005334081
Médio incompleto ou equivalente 0.005893286
Médio completo ou equivalente 0.003156167
Superior incompleto ou equivalente 0.007745619
Superior completo 0.003870122
##HISTOGRAMAS
svyhist(formula=~as.numeric(HorasSemanaisTrabalho), design=Base, main="Histograma", xlab="Número de Horas Normalmente Trabalhadas")#formula: variavel
#main,xlab:titulo, subtitulo
svyhist(formula=~as.numeric(HorasSemanaisTrabalho), design=Base, freq=TRUE, main="Histograma", xlab="Número de Horas Normalmente Trabalhadas")##BOXBLOTS
svyboxplot(formula=HorasSemanaisTrabalho~1, design=Base, main="Boxplot do Número de Horas Normalmente Trabalhadas")#~ ->variavel de agrupamento, caso queira sem, colocar 1
svyboxplot(formula=HorasSemanaisTrabalho~Sexo, design=Base, main="Boxplot do Número de Horas Normalmente Trabalhadas por Sexo")#APENAS OUTLIERS MAIS EXTREMOS
svyboxplot(formula=HorasSemanaisTrabalho~Sexo, design=Base, all.outliers=TRUE, main="Boxplot do Número de Horas Normalmente Trabalhadas por Sexo")##TESTES DE HIPOTESES
# Teste 1: Horas trabalhadas por sexo
pns_teste1 <- subset(Base, !is.na(HorasSemanaisTrabalho) & !is.na(Sexo)) #filtrando daados validos
teste1 <- svyttest(HorasSemanaisTrabalho ~ Sexo, design = pns_teste1)
teste1
Design-based t-test
data: HorasSemanaisTrabalho ~ Sexo
t = -39.962, df = 7451, p-value < 2.2e-16
alternative hypothesis: true difference in mean is not equal to 0
95 percent confidence interval:
-5.686913 -5.155071
sample estimates:
difference in mean
-5.420992
Interpretacao: Diferenca media de 5.42 horas
Conclusao: Ha diferenca significativa de horas trabalhadas entre sexos.
#REGRESSÃO LINEAR
modeloLin <- svyglm(formula=Rendimento_bruto~Nível_Instrução+Raça+Idade, design=Base)
summary(object=modeloLin)
Call:
svyglm(formula = Rendimento_bruto ~ Nível_Instrução + Raça +
Idade, design = Base)
Survey design:
subset(base, vars[[nome_coluna]] %in% c(0, 1))
Coefficients:
Estimate Std. Error
(Intercept) -658.631 93.557
Nível_InstruçãoFundamental incompleto ou equivalente 420.194 36.784
Nível_InstruçãoFundamental completo ou equivalente 856.624 44.061
Nível_InstruçãoMédio incompleto ou equivalente 1050.529 50.686
Nível_InstruçãoMédio completo ou equivalente 1311.015 55.812
Nível_InstruçãoSuperior incompleto ou equivalente 1790.372 70.420
Nível_InstruçãoSuperior completo 4077.126 106.676
RaçaPreta -670.738 46.133
RaçaAmarela -117.818 214.938
RaçaParda -608.849 41.372
RaçaIndígena -746.442 97.356
RaçaIgnorado 363.825 1324.524
Idade 39.837 1.872
t value Pr(>|t|)
(Intercept) -7.040 2.10e-12 ***
Nível_InstruçãoFundamental incompleto ou equivalente 11.423 < 2e-16 ***
Nível_InstruçãoFundamental completo ou equivalente 19.442 < 2e-16 ***
Nível_InstruçãoMédio incompleto ou equivalente 20.726 < 2e-16 ***
Nível_InstruçãoMédio completo ou equivalente 23.490 < 2e-16 ***
Nível_InstruçãoSuperior incompleto ou equivalente 25.424 < 2e-16 ***
Nível_InstruçãoSuperior completo 38.220 < 2e-16 ***
RaçaPreta -14.539 < 2e-16 ***
RaçaAmarela -0.548 0.584
RaçaParda -14.716 < 2e-16 ***
RaçaIndígena -7.667 1.98e-14 ***
RaçaIgnorado 0.275 0.784
Idade 21.285 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for gaussian family taken to be 5416107)
Number of Fisher Scoring iterations: 2
##Intervalos de confiança 95%
ic_linear<-confint(modeloLin)
coef_df <- data.frame(
Variavel = names(coef(modeloLin)),
Coeficiente = as.numeric(coef(modeloLin)),
IC_Inf = ic_linear[, 1],
IC_Sup = ic_linear[, 2]
)
coef_df Variavel
(Intercept) (Intercept)
Nível_InstruçãoFundamental incompleto ou equivalente Nível_InstruçãoFundamental incompleto ou equivalente
Nível_InstruçãoFundamental completo ou equivalente Nível_InstruçãoFundamental completo ou equivalente
Nível_InstruçãoMédio incompleto ou equivalente Nível_InstruçãoMédio incompleto ou equivalente
Nível_InstruçãoMédio completo ou equivalente Nível_InstruçãoMédio completo ou equivalente
Nível_InstruçãoSuperior incompleto ou equivalente Nível_InstruçãoSuperior incompleto ou equivalente
Nível_InstruçãoSuperior completo Nível_InstruçãoSuperior completo
RaçaPreta RaçaPreta
RaçaAmarela RaçaAmarela
RaçaParda RaçaParda
RaçaIndígena RaçaIndígena
RaçaIgnorado RaçaIgnorado
Idade Idade
Coeficiente IC_Inf
(Intercept) -658.63051 -842.02949
Nível_InstruçãoFundamental incompleto ou equivalente 420.19439 348.08821
Nível_InstruçãoFundamental completo ou equivalente 856.62422 770.25288
Nível_InstruçãoMédio incompleto ou equivalente 1050.52940 951.17111
Nível_InstruçãoMédio completo ou equivalente 1311.01518 1201.60830
Nível_InstruçãoSuperior incompleto ou equivalente 1790.37207 1652.32816
Nível_InstruçãoSuperior completo 4077.12608 3868.01043
RaçaPreta -670.73813 -761.17246
RaçaAmarela -117.81785 -539.15746
RaçaParda -608.84930 -689.95080
RaçaIndígena -746.44221 -937.28735
RaçaIgnorado 363.82458 -2232.61649
Idade 39.83676 36.16794
IC_Sup
(Intercept) -475.23153
Nível_InstruçãoFundamental incompleto ou equivalente 492.30056
Nível_InstruçãoFundamental completo ou equivalente 942.99555
Nível_InstruçãoMédio incompleto ou equivalente 1149.88770
Nível_InstruçãoMédio completo ou equivalente 1420.42206
Nível_InstruçãoSuperior incompleto ou equivalente 1928.41597
Nível_InstruçãoSuperior completo 4286.24173
RaçaPreta -580.30380
RaçaAmarela 303.52177
RaçaParda -527.74780
RaçaIndígena -555.59706
RaçaIgnorado 2960.26564
Idade 43.50557
#REGRESSÃO LOG
modeloLog <- svyglm(formula=Doença~Sexo+Raça+Idade, design=Base, family="binomial")
summary(object=modeloLog)
Call:
svyglm(formula = Doença ~ Sexo + Raça + Idade, design = Base,
family = "binomial")
Survey design:
subset(base, vars[[nome_coluna]] %in% c(0, 1))
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -3.1782216 0.0290104 -109.554 < 2e-16 ***
SexoMulher 0.3486645 0.0158377 22.015 < 2e-16 ***
RaçaPreta 0.0013958 0.0277866 0.050 0.960
RaçaAmarela 0.2138631 0.2078184 1.029 0.303
RaçaParda -0.0919759 0.0189782 -4.846 1.28e-06 ***
RaçaIndígena -0.0042507 0.1099412 -0.039 0.969
RaçaIgnorado -0.0345987 0.6965349 -0.050 0.960
Idade 0.0559756 0.0004777 117.176 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1.006224)
Number of Fisher Scoring iterations: 4
# Obter os coeficientes (log-odds) e seus intervalos de confiança
log_odds <- coef(modeloLog)
ic_log_odds <- confint(modeloLog)
# Exponenciar para obter os Odds Ratios e seus ICs, evitando valores negativos
resultado_or <- data.frame(
Variavel = names(log_odds),
OR = round(exp(log_odds), 3),
IC_Inferior = round(exp(ic_log_odds[, 1]), 3),
IC_Superior = round(exp(ic_log_odds[, 2]), 3)
)
print(resultado_or) Variavel OR IC_Inferior IC_Superior
(Intercept) (Intercept) 0.042 0.039 0.044
SexoMulher SexoMulher 1.417 1.374 1.462
RaçaPreta RaçaPreta 1.001 0.948 1.057
RaçaAmarela RaçaAmarela 1.238 0.824 1.861
RaçaParda RaçaParda 0.912 0.879 0.947
RaçaIndígena RaçaIndígena 0.996 0.803 1.235
RaçaIgnorado RaçaIgnorado 0.966 0.247 3.784
Idade Idade 1.058 1.057 1.059
##Índice de Gini
##FAZER ESTE POR ULTIMO
dadosPNS <- convey_prep(design=Base)
giniBR <- svygini(formula=~Rendimento_bruto, design=dadosPNS, na.rm=TRUE)
giniBR gini SE
Rendimento_bruto 0.50802 0.0056
Rendimento_bruto
Rendimento_bruto 0.01094457
giniUF <- svyby(formula=~Rendimento_bruto, by=~V0001, design=dadosPNS, FUN=svygini, na.rm=TRUE) #esse demora um pouco..
giniUF V0001 Rendimento_bruto se
Rondônia Rondônia 0.4359608 0.013213467
Acre Acre 0.4776399 0.011391624
Amazonas Amazonas 0.5124953 0.020344196
Roraima Roraima 0.4999963 0.017543810
Pará Pará 0.4835167 0.012171980
Amapá Amapá 0.5038926 0.013597470
Tocantins Tocantins 0.5081798 0.018013146
Maranhão Maranhão 0.5024393 0.012978459
Piauí Piauí 0.5567044 0.022037635
Ceará Ceará 0.5324006 0.023822933
Rio Grande do Norte Rio Grande do Norte 0.5260489 0.020666344
Paraíba Paraíba 0.5409785 0.019179843
Pernambuco Pernambuco 0.4934237 0.016054071
Alagoas Alagoas 0.4748417 0.017280017
Sergipe Sergipe 0.5341702 0.017287821
Bahia Bahia 0.5330190 0.017696438
Minas Gerais Minas Gerais 0.4540088 0.014635889
Espírito Santo Espírito Santo 0.4738796 0.014263516
Rio de Janeiro Rio de Janeiro 0.5154285 0.017859004
São Paulo São Paulo 0.5101548 0.013463160
Paraná Paraná 0.4535436 0.011302062
Santa Catarina Santa Catarina 0.3941068 0.009037718
Rio Grande do Sul Rio Grande do Sul 0.4683207 0.012028102
Mato Grosso do Sul Mato Grosso do Sul 0.4790102 0.019008431
Mato Grosso Mato Grosso 0.4503289 0.014585596
Goiás Goiás 0.4188121 0.009881148
Distrito Federal Distrito Federal 0.5396514 0.011241314