Os modelos de regressão com uma mistura de variáveis quantitativas e qualitativas são chamados de modelos de análise de covariância (ANCOVA). Eles são uma extensão dos modelos ANOVA no sentido de que fornecem um método de controle estatístico dos efeitos das variáveis independentes quantitativas, chamadas covariáveis ou variáveis de controle, em um modelo que inclui tanto variáveis quantitativas quanto qualitativas ou binários.
Vamos seguir nosso exemplo prático anterior com base no salário médio (em dólares) de professores de escolas públicas nos EUA para o ano escolar de 2005–2006. Não obstante, agora vamos incluir ao modelo a variável gastos das autoridades locais com a escola pública.
Repassando alguns passos.
Primeiro, importamos o conjunto de dados.
library(readxl)
exercicio_anova <- read_excel("D:/Documentos/Fabio/OneDrive/Área de Trabalho/ADRI_exercício/exercicio_anova.xlsx", col_types = c("text", "text", "numeric", "numeric"))
print(exercicio_anova)## # A tibble: 51 × 4
## Estado Regiao Salario Gastos
## <chr> <chr> <dbl> <dbl>
## 1 Connecticut Nordeste 60822 12436
## 2 Illinois Centro-Norte 58246 9275
## 3 Indiana Centro-Norte 47831 8935
## 4 Iowa Centro-Norte 43130 7807
## 5 Kansas Centro-Norte 43334 8373
## 6 Maine Nordeste 41596 11285
## 7 Massachusetts Nordeste 58624 12596
## 8 Michigan Centro-Norte 54895 9880
## 9 Minnesota Centro-Norte 49634 9675
## 10 Missouri Centro-Norte 41839 7840
## # ℹ 41 more rowsComo se pode observar, o conjunto de dados apresenta 51 observações classificadas em três regiões geográficas:
Para realizar a análise de variancia é preciso definir as variáveis binárias para representar cada região.
Para definir as variáveis binárias que representam cada região, é necessário selecionar uma delas como região de referência.
Para o exercício, considere a seguinte classificação:
Para criar as variávies binárias no R, execute o seguinte código:
library(dplyr)
exercicio_anova <- exercicio_anova %>%
mutate(D2 = case_when(
Regiao == "Nordeste" | Regiao == "Centro-Norte" ~ 1,
TRUE ~ 0
))
exercicio_anova <- exercicio_anova %>%
mutate(D3 = case_when(
Regiao == "Sul" ~ 1,
TRUE ~ 0
))
print(exercicio_anova)## # A tibble: 51 × 6
## Estado Regiao Salario Gastos D2 D3
## <chr> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 Connecticut Nordeste 60822 12436 1 0
## 2 Illinois Centro-Norte 58246 9275 1 0
## 3 Indiana Centro-Norte 47831 8935 1 0
## 4 Iowa Centro-Norte 43130 7807 1 0
## 5 Kansas Centro-Norte 43334 8373 1 0
## 6 Maine Nordeste 41596 11285 1 0
## 7 Massachusetts Nordeste 58624 12596 1 0
## 8 Michigan Centro-Norte 54895 9880 1 0
## 9 Minnesota Centro-Norte 49634 9675 1 0
## 10 Missouri Centro-Norte 41839 7840 1 0
## # ℹ 41 more rowsObserve que duas novas colunas, D2 e D3, foram adicionadas ao conjunto de dados. Quando ambas assumem o valor zero, a observação corresponde à região de referência, representada por D1.
Para que o R interprete uma variável numérica como qualitativa, é importante convertê-la em um fator. Essa transformação é necessária porque a ANOVA compara categorias, e não valores contínuos. Portanto, o modelo precisa reconhecer que os valores representam grupos distintos, e não uma escala numérica.
exercicio_anova$D2 <- as.factor(exercicio_anova$D2)
exercicio_anova$D3 <- as.factor(exercicio_anova$D3)
print(exercicio_anova)## # A tibble: 51 × 6
## Estado Regiao Salario Gastos D2 D3
## <chr> <chr> <dbl> <dbl> <fct> <fct>
## 1 Connecticut Nordeste 60822 12436 1 0
## 2 Illinois Centro-Norte 58246 9275 1 0
## 3 Indiana Centro-Norte 47831 8935 1 0
## 4 Iowa Centro-Norte 43130 7807 1 0
## 5 Kansas Centro-Norte 43334 8373 1 0
## 6 Maine Nordeste 41596 11285 1 0
## 7 Massachusetts Nordeste 58624 12596 1 0
## 8 Michigan Centro-Norte 54895 9880 1 0
## 9 Minnesota Centro-Norte 49634 9675 1 0
## 10 Missouri Centro-Norte 41839 7840 1 0
## # ℹ 41 more rowsAgora o conjunto de dados está pronto para rodar o modelo de regressão.
Suponhamos que desejamos verificar se o salário médio anual dos professores de escolas públicas varia conforme as diferentes regiões em que o país foi dividido e os gastos das autoridades locais com a escola pública.
Para tanto, podemos aplicar a ANCOVA para comparar dois ou mais valores médios com base em variáveis qualitativas e quantitativas.
Vamos admitir o seguinte modelo:
\[Y_{i} = \beta_{1} + \beta_{2}D_{2i} + \beta_{3}D_{3i} + \beta_{4}X_{i} + u_{i}\]
onde,
\(Y_{i}\) = salário médio dos professores da rede pública no Estado i;
\(D_{2i}\) = 1 se o estado pertencer a região Nordeste ou Norte-central; = 0 se estiver situado em outras regioões do país;
\(D_{3i}\) = 1 se o estado pertencer a região Sul; 0 se estiver localizado em outras regiões;
\(X_{i}\) = gastos com escolas públicas por aluno;
\(u_{i}\) = termo do erro.
modelo_anova <- lm(Salario ~ D2 + D3 + Gastos, data = exercicio_anova)
summary(modelo_anova)##
## Call:
## lm(formula = Salario ~ D2 + D3 + Gastos, data = exercicio_anova)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10556 -2471 106 2066 15084
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 28694.9180 3262.5213 8.795 1.70e-11 ***
## D21 -2954.1268 1862.5756 -1.586 0.1194
## D31 -3112.1948 1819.8725 -1.710 0.0938 .
## Gastos 2.3404 0.3592 6.515 4.45e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4905 on 47 degrees of freedom
## Multiple R-squared: 0.4977, Adjusted R-squared: 0.4656
## F-statistic: 15.52 on 3 and 47 DF, p-value: 3.762e-07O modelo ANCOVA sugere que os gastos com educação pública por aluno têm um impacto forte e positivo nos salários do professores, enquanto as diferenças salariais entre regiões não são estatisticamente significativas. Isso pode indicar que o efeito regional sobre os salários precisa ser analisado com mais profundidade ou incluir outras variáveis que possam influenciar essas diferenças.
Vamos analisar agora um outro exemplo ANCOVA, levando em consideração a renda, a idade, o genêro e o grau de escolaridade.
Primeiro, importamos o conjunto de dados que apresenta informações referentes ao rendimento semanal de trabalhadores indianos de uma cidade industrial em 1990.
library(readxl)
exercicio_india <- read_excel("D:/Documentos/Fabio/OneDrive/Área de Trabalho/ADRI_exercício/exercicio_india.xlsx")
print(exercicio_india)## # A tibble: 114 × 7
## RS Idade DE2 DE3 DE4 DPT Dgen
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 120 57 0 0 0 0 0
## 2 224 48 0 0 1 1 0
## 3 132 38 0 0 0 0 0
## 4 75 27 0 1 0 0 0
## 5 111 23 0 1 0 0 1
## 6 127 22 0 1 0 0 0
## 7 30 18 0 0 0 0 0
## 8 24 12 0 0 0 0 0
## 9 119 38 0 0 0 1 0
## 10 75 55 0 0 0 0 0
## # ℹ 104 more rowsOnde,
\(RS\) = renda semanal em rúpias \(Idade\) = idade em anos \(D_{gen}\) = 0 para homens, 1 para mulheres; \(DPT\) = 0 para trabalhadores temporários, 1 para trabalhadores contratados com tempo indeterminado; \(DE2\) = 1 para trabalhadores com fundamental grau completo, 0 para os demais; \(DE3\) = 1 para trabalhadores com médio grau completo, 0 para o demais; \(DE4\) = 1 para trabalhadores com ensino superior completo, 0 para os demais.
exercicio_india$DE2 <- as.factor(exercicio_india$DE2)
exercicio_india$DE3 <- as.factor(exercicio_india$DE3)
exercicio_india$DE4 <- as.factor(exercicio_india$DE4)
exercicio_india$DPT <- as.factor(exercicio_india$DPT)
exercicio_india$Dgen <- as.factor(exercicio_india$Dgen)
print(exercicio_india)## # A tibble: 114 × 7
## RS Idade DE2 DE3 DE4 DPT Dgen
## <dbl> <dbl> <fct> <fct> <fct> <fct> <fct>
## 1 120 57 0 0 0 0 0
## 2 224 48 0 0 1 1 0
## 3 132 38 0 0 0 0 0
## 4 75 27 0 1 0 0 0
## 5 111 23 0 1 0 0 1
## 6 127 22 0 1 0 0 0
## 7 30 18 0 0 0 0 0
## 8 24 12 0 0 0 0 0
## 9 119 38 0 0 0 1 0
## 10 75 55 0 0 0 0 0
## # ℹ 104 more rowsPara especificar o modelo, vamos admitir como categoria de referência os trabalhadores do gênero masculino com primeiro grau incompleto e trabalho temporário.
Portanto, que satisfaz a seguinte condição:
\(D_{gen} = 0\)
\(DPT = 0\)
\(DE2 = 0\)
\(DE3 = 0\)
\(DE4 = 0\)
Dito isso, nosso objetivo é investigar de que forma os salários semanais se relacionam com a idade, o gênero, o nível de escolaridade e o tempo de serviço no emprego.
Para tanto, estimamos o seguinte modelo de regressão:
\[ln(RS_{i} = \beta_{1} + \beta_{2}Idade + \beta_{3}D_{gen} + \beta_{4}DE2 + \beta_{5}DE3 + \beta_{6}DE4 + \beta_{7}DPT + u_{i})\]
Seguindo a literatura da Economia do Trabalho, expressamos o logaritmo natural dos salários como uma função das variáveis explanatórias. Essa abordagem é justificada pelo fato de que, conforme observado, a distribuição de variáveis como os salários tende a ser assimétrica. A transformação logarítmica contribui para reduzir essa assimetria, além de mitigar problemas de heterocedasticidade nos resíduos do modelo.
Com isso tempos:
modelo_india <- lm(log(RS) ~ Idade + Dgen + DE2 + DE3 + DE4 + DPT, data = exercicio_india)
summary(modelo_india)##
## Call:
## lm(formula = log(RS) ~ Idade + Dgen + DE2 + DE3 + DE4 + DPT,
## data = exercicio_india)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.76975 -0.36798 0.04408 0.33801 1.97701
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.622107 0.169943 21.314 < 2e-16 ***
## Idade 0.029309 0.004708 6.226 9.58e-09 ***
## Dgen1 -0.666497 0.143860 -4.633 1.02e-05 ***
## DE21 0.157655 0.165072 0.955 0.34170
## DE31 0.505008 0.165203 3.057 0.00282 **
## DE41 0.623553 0.266812 2.337 0.02130 *
## DPT1 0.319470 0.126350 2.528 0.01292 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6072 on 107 degrees of freedom
## Multiple R-squared: 0.4992, Adjusted R-squared: 0.4711
## F-statistic: 17.78 on 6 and 107 DF, p-value: 3.349e-14Resultado
Considerando que o modelo satisfaça todas as premiças do MQO, os resultados mostram que o logaritmo de salários está positivamente relacionado à idade, educação e permanência no emprego, mas negativamente relacionado ao gênero.
Embora pareça não haver diferença prática nos salários semanais de trabalhadores com graus de escolaridade primário ou menor, os salários semanais são mais altos para trabalhadores com segundo grau e muito mais altos para trabalhadores com educação superior.
Lembre-se: os coeficientes das variáveis binárias devem ser interpretados como valores diferenciais da categoria de referência.
O coeficiente da variável DPT sugere que aqueles trabalhadores com contratos por tempo indeterminado ganham, em média, mais dinheiro que aqueles cujos trabalhos são temporários.
OBS: em um modelo log-linear (variáveis dependentes em forma de logaritmo e variáveis explanatórias em forma linear), o coeficiente angular de uma variável explanatória representa a semielasticidade, ou seja, fornece a variação percentual ou relativa na variável dependente para uma variação de unidade no valor da variável explanatória.
Não obstante, quando a variável explanatória é uma variável dummy, temos de ser cuidadosos. Aqui temos de tomar o antilogaritmo do coeficiente binário estimado, subtrair 1 dele e multiplicar o resultado por 100. Ou seja, para descobrirmos a variação percentual em salários semanais para aqueles trabalhadores que têm empregos por tempo indeterminado versus aqueles que têm empregos temporários, tomamos o antilogaritmo do coeficiente DPT, subtraímos 1 e então multiplicamos a diferença por 100.
Para o nosso exemplo:
(exp(modelo_india$coefficients[["DPT1"]]) - 1)*100## [1] 37.63978Para o nosso exemplo, mantendo-se constantes as demais variáveis, trabalhadores que possuem contrato por tempo indeterminado possuem, em média, uma renda semanal 37.6% maior em relação aos trabalhadores com contrato temporário.
Além disso, o estudo indica que gênero e escolaridade também têm efeitos diferencias nos ganhos semanais. Vejamos os efeitos diferencias tomados individualmente.
(exp(modelo_india$coefficients[["Dgen1"]]) - 1)*100## [1] -48.64957Portanto, mantendo-se constantes as demais variáveis, as mulheres indianas recebem, em média, 48,64% a menos por semana em comparação aos homens.
(exp(modelo_india$coefficients[["DE31"]]) - 1)*100## [1] 65.69983Mantendo-se constantes as demais variáveis, os trabalhadores indianos com ensino médio completo recebem, em média, 65,7% a mais por semana em comparação aos demais níveis de escolaridade.
(exp(modelo_india$coefficients[["DE41"]]) - 1)*100## [1] 86.5545Já, mantendo-se constantes as demais variáveis, os trabalhadores indianos com ensino superior completo recebem, em média, 86,55% a mais por semana em comparação aos demais níveis de escolaridade.
modelo_india2 <- lm(log(RS) ~ Idade + Dgen + DE2 + DE3 + DE4 + Dgen*DE2 + Dgen*DE3 + Dgen*DE4 + DPT, data = exercicio_india)
summary(modelo_india2)##
## Call:
## lm(formula = log(RS) ~ Idade + Dgen + DE2 + DE3 + DE4 + Dgen *
## DE2 + Dgen * DE3 + Dgen * DE4 + DPT, data = exercicio_india)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.82026 -0.37652 0.04838 0.30159 1.95085
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.642675 0.170947 21.309 < 2e-16 ***
## Idade 0.029575 0.004735 6.247 9.31e-09 ***
## Dgen1 -0.804700 0.167351 -4.808 5.15e-06 ***
## DE21 0.064557 0.177382 0.364 0.7166
## DE31 0.432772 0.177745 2.435 0.0166 *
## DE41 0.449705 0.364887 1.232 0.2206
## DPT1 0.332894 0.127692 2.607 0.0105 *
## Dgen1:DE21 0.655485 0.489911 1.338 0.1838
## Dgen1:DE31 0.461928 0.486463 0.950 0.3445
## Dgen1:DE41 0.402677 0.525642 0.766 0.4454
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6079 on 104 degrees of freedom
## Multiple R-squared: 0.5121, Adjusted R-squared: 0.4698
## F-statistic: 12.13 on 9 and 104 DF, p-value: 6.446e-13