library(openintro)
library(tidyverse)
library(ggbeeswarm)
library(modelr)
library(broom)
theme_set(theme_bw())
O exercício é sobre regressão com múltiplas variáveis, nele temos um conjunto de dados sobre professores e desejamos encontrar o modelo linear que melhor represente a nossa variável resposta em função de todas as outras.
Os dados que usaremos estão detalhados nessa página do OpenIntro.
profs = read_csv("https://www.openintro.org/stat/data/evals.csv") %>%
select(score, age, gender, bty_avg, pic_outfit, pic_color)
## Parsed with column specification:
## cols(
## .default = col_double(),
## rank = col_character(),
## ethnicity = col_character(),
## gender = col_character(),
## language = col_character(),
## cls_level = col_character(),
## cls_profs = col_character(),
## cls_credits = col_character(),
## pic_outfit = col_character(),
## pic_color = col_character()
## )
## See spec(...) for full column specifications.
glimpse(profs)
## Rows: 463
## Columns: 6
## $ score <dbl> 4.7, 4.1, 3.9, 4.8, 4.6, 4.3, 2.8, 4.1, 3.4, 4.5, 3.8, 4.5…
## $ age <dbl> 36, 36, 36, 36, 59, 59, 59, 51, 51, 40, 40, 40, 40, 40, 40…
## $ gender <chr> "female", "female", "female", "female", "male", "male", "m…
## $ bty_avg <dbl> 5.000, 5.000, 5.000, 5.000, 3.000, 3.000, 3.000, 3.333, 3.…
## $ pic_outfit <chr> "not formal", "not formal", "not formal", "not formal", "n…
## $ pic_color <chr> "color", "color", "color", "color", "color", "color", "col…
Acima podemos ter uma visão geral dos nossos dados observando as variáveis disponíveis para análise. As nossas variáveis são score que representa a pontuação de um professor, age que é a idade dos professores, gender que é o gênero do professor, bty_avg que é uma pontuação quanto a beleza do professor, pic_outfit que indica se a foto do professor é formal ou não e pic_color que indica se a foto é colorida ou não.
Abaixo poderemos ver algumas visualizações gráficas que relacionam o score com os outros atributos do professor e com isso poderemos ter uma idea inicial das relações entre os atributos.
profs %>%
ggplot(aes(x = age, y = score)) +
geom_count()
profs %>%
ggplot(aes(x = bty_avg, y = score)) +
geom_count()
profs %>%
ggplot(aes(x = pic_outfit, y = score)) +
geom_count()
profs %>%
ggplot(aes(x = pic_color, y = score)) +
geom_count()
Os gráficos que comparam score com a idade e com o nível de beleza respectivamente, não deixam muito claro como o score varia em função de ambas as variáveis, mas quanto ao gráfico relativo a idade, podemos destacar que a concentração de professores é maior abaixo de 60 anos. Os gráficos relativos a formalidade da foto e a cor da foto respectivamente tem resultados bem semelhantes, no gráfico que analisa a formalidade a amplitude dos scores quando a foto é informal é maior, assim como a concentração de fotos informais é maior, enquanto que no gráfico relativo a cor da foto a amplitude dos scores quando a foto é colorida é maior e a concentração de fotos coloridas também é maior.
Abaixo vamos visualizar o gráfico de densidade da variável score que na nossa regressão será nossa variável de respostas. Com a análise do gráfico abaixo poderemos ter uma ideia de por volta de quais valores se concentram os scores dos professores.
profs %>%
ggplot(aes(x = score)) +
geom_density(fill = "darkred") +
labs(title = "Concentração dos scores", x = "Score", y = "Densidade")
Podemos ver no gráfico acima que os scores dos professores se concentram mais em torno de 4.5, que é uma nota bem alta, mais adiante vamos ver que fatores influênciam no resultado do score de um determinado professor.
Abaixo podemos ver um sumário estatístico com a média, mediana e desvio padrão da variável na qual estamos interessados: score.
profs %>%
summarise(média_scores = mean(score), mediana_score = median(score),
desvio_padrao = sd(score))
## # A tibble: 1 x 3
## média_scores mediana_score desvio_padrao
## <dbl> <dbl> <dbl>
## 1 4.17 4.3 0.544
Com base nos resultados acima, podemos dizer que não há muitos dados destoantes do comportamento geral da variável, a média e a mediana são muito próximas e o desvio padrão é baixo, ou seja, há pouca chance de haver pontos extremos. Além disso, podemos reforçar o que foi observado no gráfico de densidade, os scores dos professores realmente tendem a ser valores altos, se concentrando bastante acima de 4.
Queremos descobrir o efeito da variação das variáveis de preditoras (age, gender, bty_avg, pic_color, pic_outfit) sobre a variável resposta que é score.
modelo = lm(score ~ age + gender + bty_avg + pic_color + pic_outfit, profs)
tidy(modelo, conf.int = TRUE) %>% select(term, estimate, conf.low, conf.high)
## # A tibble: 6 x 4
## term estimate conf.low conf.high
## <chr> <dbl> <dbl> <dbl>
## 1 (Intercept) 4.37 3.96 4.78
## 2 age -0.00712 -0.0126 -0.00167
## 3 gendermale 0.223 0.121 0.326
## 4 bty_avg 0.0473 0.0126 0.0821
## 5 pic_colorcolor -0.212 -0.348 -0.0765
## 6 pic_outfitnot formal -0.0152 -0.148 0.118
Com os resultados calculado acima podemos ver graficamente os intervalos onde temos confiança que os coeficientes calculados para cada uma das variáveis preditivas está.
modelo %>%
tidy(conf.int = TRUE) %>%
filter(term != "(Intercept)") %>%
ggplot(aes(estimate, reorder(term, estimate),
xmin = conf.low,
xmax = conf.high,
height = 0.1)) +
geom_vline(xintercept = 0, color = "grey") +
# geom_point() +
geom_errorbarh()
Com base base no resultado podemos montar o seguinte modelo: \(score = 4,4 + 0,22.gender + 0,05.btyavg - 0,007.age - 0,21.piccolor - 0,015.picoutfit\). Com base nesse resultado podemos tirar alguma conclusões sobre as variáveis preditivas:
A variável com maior efeito sobre o score provavelmente é gender, se o professor for homem o efeito é positivo e está no intervalo entre 0.12 e 0.32.
A variável pic_color também possui um efeito quase tão grande quanto o de gender, estando no intervalo entre -0.35 e -0.07 quando a foto do professor é colorida.
A influência da age no nosso modelo é mínima e negativa, ou seja, quanto mais velho, pior tende a ser a avaliação, a influência é irrelevante quando o professor é jovem, estando em um intervalo entre -0.012 e -0.0016, mas quando ele é mais velho, essa influência pode chegar a ser de 0.315.
A influência do pic_outfit quando a foto não é formal, se existir, pode ser negativa ou positiva na avaliação, o intervalo é entre -0.15 e 0.12, ou seja, a influência sobre o resultado do modelo é muito baixa.
O intervalo do bty_avg é entre 0.012 e 0.082, ou seja, pode ser o fator que menos influencia o modelo.
Podemos reforçar a ideia de que o valor do score tende a ser altos pelo intercept que vale 4.4 quando todas as outras variáveis são 0.
glance(modelo)
## # A tibble: 1 x 11
## r.squared adj.r.squared sigma statistic p.value df logLik AIC BIC
## <dbl> <dbl> <dbl> <dbl> <dbl> <int> <dbl> <dbl> <dbl>
## 1 0.0878 0.0778 0.522 8.79 5.60e-8 6 -353. 720. 749.
## # … with 2 more variables: deviance <dbl>, df.residual <int>
Como podemos ver no resultado do \(R^2\) o nosso modelo representa muito pouco da relação, sendo \(R^2 = 0,087\), ou seja, representando cerca de 8,7% da relação que analisamos.