Correlação / Regressão
Estatística Aplicada - UFPR - Maury Melo
20 de junho de 2019
PRIMEIRA PARTE
Correlação - Pesquisa Geral
Primeiro Exercício
Um investigador deseja estudar a possível relação entre o salário (em anos completos) e o tempo de experiência (em mil reais) no cargo de gerente de agências bancárias de uma grande empresa. Os dados coletados são mostrados abaixo:
Solução: Inicialmente, os dados devem ser organizados como objetos de dados R, nesse caso como um data frame (planilha). Para isso, é necessário que a tabela acima se encontre numa estrutura tabular, na qual as colunas representam as variáveis e as linhas representam os indivíduos.
# Salario Experiencia
# 1.9307 0
# 3.1769 17
# 2.2769 8
# 3.1307 15
# 2.7769 9
# 3.0923 15
# 2.6508 8
# 2.2230 5
# 2.8538 13
# 3.2307 20
# 2.8230 11
# 1.9076 1
# 2.5384 6
# 2.5692 7
# 4.2230 23
# 4.0923 20
# 3.6000 18
# 4.7076 27
# 3.1461 11
# 2.9923 10
# 4.7461 29
# 4.1153 23
# 2.3615 4
# 4.0923 22
# 4.5076 25
# 2.9076 9
# 4.4846 25
Salario=c(1.9307,3.1769,2.2769,3.1307,2.7769,3.0923,2.6508,2.2230,2.8538,3.2307,2.8230,1.9076,2.5384,2.5692,4.2230,4.0923,3.6000,4.7076,3.1461,2.9923,4.7461,4.1153,2.3615,4.0923,4.5076,2.9076,4.4846)
Experiencia=c(0,17,8,15,9,15,8,5,13,20,11,1,6,7,23,20,18,27,11,10,29,23,4,22,25,9,25)
cbind(Salario,Experiencia)## Salario Experiencia
## [1,] 1.9307 0
## [2,] 3.1769 17
## [3,] 2.2769 8
## [4,] 3.1307 15
## [5,] 2.7769 9
## [6,] 3.0923 15
## [7,] 2.6508 8
## [8,] 2.2230 5
## [9,] 2.8538 13
## [10,] 3.2307 20
## [11,] 2.8230 11
## [12,] 1.9076 1
## [13,] 2.5384 6
## [14,] 2.5692 7
## [15,] 4.2230 23
## [16,] 4.0923 20
## [17,] 3.6000 18
## [18,] 4.7076 27
## [19,] 3.1461 11
## [20,] 2.9923 10
## [21,] 4.7461 29
## [22,] 4.1153 23
## [23,] 2.3615 4
## [24,] 4.0923 22
## [25,] 4.5076 25
## [26,] 2.9076 9
## [27,] 4.4846 25Note que são considerados 27 pares de observações correspondentes à variável resposta Salário e à variável explicativa Experiência, para cada um dos gerentes da empresa.
Análise Exploratória
summary(Salario)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1.908 2.610 3.092 3.228 4.092 4.746var(Salario)## [1] 0.7368296sd(Salario)## [1] 0.8583878Uma maneira fácil de obter algumas estatísticas descritivas das variáveis em estudo é através do comando summary(), que retorna as estatísticas mínimo, quartis, média e máximo. Para medir a variabilidade, foram utilizadas as funções var() e sd() para obter a variância e o desvio padrão.
Diagrama de Dispersão
Para verificar a existência de alguma relação entre Salário e Experiência, será construído um Diagrama de Dispersão para as duas variáveis:
m = lm(Experiencia~Salario)
plot(Experiencia~Salario)
abline(m)
O gráfico de dispersão é usado para verificar se existe relação de causa e efeito entre duas variáveis de natureza quantitativa (variáveis que podem ser medidas ou contadas). Isto não prova que uma variável afeta outra variável, mas determina se existe relação e qual a intensidade da relação entre elas.No nosso caso pode-se observar que existe uma correlação linear, bem como identificar que o valor da variável “Salário” aumenta com o aumento do valor da variável “experiência” - neste caso, diz–se que as variáveis são positivamente correlacionadas.
O maior objetivo do boxplot é verificar a distribuição dos dados.A dispersão é representada pela amplitude do gráfico, que pode ser calculada como máximo valor – mínimo valor. Quanto maior for a amplitude, maior a variação nos dados.O retângulo contém 50% dos valores do conjunto de dados. A posição da linha mediana no retângulo informa sobre a assimetria da distribuição.Veja a seguir o boxplot de nosso exemplo.
O Boxplot fornece informação sobre as seguintes características do conjunto de dados: localização, dispersão, assimetria, comprimento da cauda e outliers (medidas discrepantes). O centro da distribuição é indicado pela linha da mediana, no centro da caixa. A dispersão é representada pela amplitude do gráfico, que pode ser calculada como máximo valor – mínimo valor. Quanto maior for a amplitude, maior a variação nos dados.
par(mfrow=c(1, 2)) # divide graph area in 2 columns
boxplot(Experiencia,data = InsectSprays, col = "lightgray", main="Experiência")
boxplot(Salario, data = InsectSprays, col = "lightgray", main="Salário")
Uma distribuição simétrica teria a mediana no centro do retângulo. Se a mediana é próxima de Q1, então, os dados são positivamente assimétricos. Se a mediana é próxima de Q3 os dados são negativamente assimétricos.Os outliers em um box plot aparecem como pontos ou asteriscos fora das “linhas” desenhada (em nosso caso não existe outlier).Uma utilidade muito importante do Boxplot é na comparação gráfica de dois ou mais grupos. Nesse caso, o Boxplot é preferível ao histograma. Essa comparação pode ser feita desenhando-se os Boxplots para cada conjunto de dados, paralelamente, em um mesmo gráfico. Observe no gráfico acima que há indicativos de que a variância dos erros é constante, e não há evidências de que os erros não sigam a distribuição Normal e percebe-se de que os erros não são correlacionados entre sí.
Para efeito de exemplificar um caso de boxplot com outlier, segue exemplo abaixo:
Histograma
Um histograma é uma representação visual da distribuição de um conjunto de dados. Como tal, a forma de um histograma é sua característica mais óbvia e informativa: ele permite que você veja facilmente onde uma quantidade relativamente grande dos dados está situada e onde há muito poucos dados a serem encontrados (Verzani 2004). Em outras palavras, você pode ver onde o meio está em sua distribuição de dados, quão próximos estão os dados em torno deste meio e, quando possível, os outliers podem ser encontrados. Exatamente por tudo isso, os histogramas são uma ótima maneira de conhecer seus dados. Este tipo de gráfico adequado para utilização com dados numéricos e se busca aprender sobre sua distribuição. No caso de densidade, a frequência relativa do intervalo i, (fri), é representada pela área de um retângulo que é colocado acima do ponto médio da classe i. Consequentemente, a área total do histograma (igual a soma das áreas de todos os retângulos) será igual a 1. Assim, ao construir o histograma, cada retângulo deverá ter área proporcional à frequência relativa (ou à frequência absoluta, o que é indiferente) correspondente.
par(mfrow=c(1, 2)) # divide graph area in 2 columns
hist(Salario,col="orange", border="black", prob=TRUE, xlab="MPG", ylab="Densidade", main="Salário")
abline(v = mean(Salario), col = "blue")
hist(Experiencia,col="green", border="black", prob=TRUE, xlab="MPG", ylab="Densidade", main="Experiência")
abline(v = mean(Experiencia), col = "red")
No histograma olhamos informações de localização e de variação dos dados. Por exemplo, qual o centro (a média), qual o mínimo e o máximo, qual a amplitude (variação) - O objetivo é conhecer um pouco sobre os dados que estamos trabalhando e também identificar a possível presença de outliers. Por essas razões, o histograma é uma das primeiras ferramentas a ser utilizado na chamada análise descritiva ou exploratória de dados.
Por meio do gráfico podemos observar o comportamento dos pesos e examinar características como: simetria da distribuição (nesse caso podemos entender como distribuição simétrica), região onde há maior concentração de valores, o centro da distribuição (conforme reta vertical), a dispersão dos valores ao redor de uma medida central. Também é possível perceber que não existe nenhum outlier .
Correlação Linear de Person
Os coeficientes de correlação são usados nas estatísticas para medir o quão forte é uma relação entre duas variáveis. Existem vários tipos de coeficientes de correlação: a correlação de Pearson (também chamada de R de Pearson) é um coeficiente de correlação comumente usado na regressão linear.
O coeficiente de correlação de Pearson é uma medida da força da relação linear entre duas variáveis. É referido como correlação de Pearson ou simplesmente como o coeficiente de correlação. Se a relação entre as variáveis não for linear, então o coeficiente de correlação não representa adequadamente a força da relação entre as variáveis. função \(cor.test()\) pode ser usada para calcular o nível de significância da correlação. Ele testa a associação entre amostras pareadas usando os métodos pearson, kendall ou spearman.
Basicamente, uma correlação momento-produto Pearson tenta traçar uma linha de melhor ajuste através dos dados de duas variáveis, e o coeficiente de correlação de Pearson, r, indica a que distância todos esses pontos de dados estão para essa linha de melhor ajuste bem os pontos de dados se encaixam nesse novo modelo / linha de melhor ajuste).
Um coeficiente de correlação de 1 significa que, para cada aumento positivo em uma variável, há um aumento positivo de uma proporção fixa na outra. Por exemplo, os tamanhos dos sapatos aumentam em correlação (quase) perfeita com o comprimento do pé.Um coeficiente de correlação de -1 significa que, para cada aumento positivo em uma variável, há uma diminuição negativa de uma proporção fixa na outra. Por exemplo, a quantidade de gás em um tanque diminui em correlação (quase) perfeita com a velocidade.Zero significa que, para cada aumento, não há um aumento positivo ou negativo. Os dois simplesmente não estão relacionados.
O coeficiente de correlação de Pearson, r, pode ter um intervalo de valores de +1 a -1. Um valor de 0 indica que não há associação entre as duas variáveis. Um valor maior que 0 indica uma associação positiva; isto é, à medida que o valor de uma variável aumenta, o mesmo acontece com o valor da outra variável. Um valor menor que 0 indica uma associação negativa; isto é, à medida que o valor de uma variável aumenta, o valor da outra variável diminui. Isso é mostrado na figura abaixo:
Quanto mais forte a associação das duas variáveis, mais próximo o coeficiente de correlação de Pearson, r, será para +1 ou -1, dependendo se a relação é positiva ou negativa, respectivamente. Alcançar um valor de +1 ou -1 significa que todos os seus pontos de dados estão incluídos na linha de melhor ajuste - não há pontos de dados que mostrem qualquer variação fora dessa linha. Valores para r entre +1 e -1 (por exemplo, r = 0,8 ou -0,4) indicam que há variação ao redor da linha de melhor ajuste. Quanto mais próximo o valor de r a 0, maior a variação em torno da linha de melhor ajuste. Relacionamentos diferentes e seus coeficientes de correlação são mostrados no diagrama abaixo:
options(scipen = 999)
cor.test(Experiencia,Salario)##
## Pearson's product-moment correlation
##
## data: Experiencia and Salario
## t = 20.1, df = 25, p-value < 0.00000000000000022
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.9353430 0.9866043
## sample estimates:
## cor
## 0.9704256Como o Valor P do teste (p-value < 2.2e-16) é bem pequeno, conclui-se que o valor do Coeficiente de Correlação Linear de Pearson tem significância Estatística.
O símbolo da correlação de Pearson é “ρ” quando é medido na população e “r” quando é medido em uma amostra. Como estaremos lidando quase exclusivamente com amostras, usaremos r para representar a correlação de Pearson, a menos que indicado de outra forma.O r de Pearson pode variar de -1 a 1. Um r de -1 indica um relacionamento linear negativo perfeito entre variáveis, um r de 0 indica nenhuma relação linear entre variáveis, e um r de 1 indica um perfeito relacionamento linear positivo entre variáveis. A figura abaixo mostra um gráfico de dispersão para o qual r = 1 e para r = -1.
Para testar se suas duas variáveis formam um relacionamento linear, basta plotar um gráfico de dispersão e inspecionar visualmente a forma do gráfico. No diagrama abaixo, você encontrará alguns exemplos diferentes de um relacionamento linear e alguns relacionamentos não lineares. Não é apropriado analisar uma relação não linear usando uma correlação de Pearson. A correlação de Pearson determina o grau em que um relacionamento é linear. Em outras palavras, determina se existe um componente linear de associação entre duas variáveis contínuas.
Teste da Normalidade de Shapiro
Os testes de normalidade são complementares à avaliação gráfica da normalidade.O teste de Shapiro compara os escores da amostra a um conjunto de escores normalmente distribuídos com a mesma média e desvio padrão; a hipótese nula é que “a distribuição da amostra é normal”. Se o teste for significativo, a distribuição não é normal. Para amostras pequenas, os testes de normalidade têm pouco poder para rejeitar a hipótese nula e, portanto, amostras pequenas geralmente passam em testes de normalidade. Para amostras grandes, resultados significativos seriam obtidos mesmo no caso de um pequeno desvio da normalidade, embora esse pequeno desvio não afete os resultados de um teste paramétrico. O teste de Shapiro-Wilk é baseado na correlação entre os dados e os escores normais correspondentes e fornece uma boa opção para testar a normalidade dos dados. O poder é a medida mais frequente do valor de um teste de normalidade - a capacidade de detectar se uma amostra provém de uma distribuição não normal.
shapiro.test(Salario)##
## Shapiro-Wilk normality test
##
## data: Salario
## W = 0.93433, p-value = 0.08831shapiro.test(Experiencia)##
## Shapiro-Wilk normality test
##
## data: Experiencia
## W = 0.96183, p-value = 0.4064
A Hipótese nula do teste de Shapiro é que a população é normalmente distribuída. Assim, por um lado, se o valor p for menor que o nível alfa escolhido, a hipótese nula é rejeitada e há evidências de que os dados testados não são normalmente distribuídos. Por outro lado, se o valor p for maior que o nível alfa escolhido, a hipótese nula de que os dados vieram de uma população normalmente distribuída não pode ser rejeitada (por exemplo, para um nível alfa de 0,05 - um conjunto de dados com valor inferior a 0,05 rejeita a hipótese nula de que os dados são de uma população normalmente distribuída). Como a maioria dos testes de significância estatística, se o tamanho da amostra for suficientemente grande, este teste pode detectar desvios mesmo triviais da hipótese nula (isto é, embora possa haver algum efeito estatisticamente significativo, pode ser muito pequeno para ter qualquer significado prático); portanto, uma investigação adicional do tamanho do efeito é normalmente aconselhável.
SEGUNDA PARTE
Análise de Regressão
Análise de regressão é uma técnica estatística utilizada para investigar a relação existente entre variáveis através da construção de uma equação (um modelo).Iniciemos com o mesmo exemplo anterior, onde um investigador deseja estudar a possível relação entre o salário (em anos completos) e o tempo de experiência (em mil reais) no cargo de gerente de agências bancárias de uma grande empresa.
Inicialmente, os dados devem ser organizados como objetos de dados R, nesse caso como um data frame (planilha). Para isso, é necessário que a tabela acima se encontre numa estrutura tabular, na qual as colunas representam as variáveis e as linhas representam os indivíduos - como feito anteriormente na Análise de Correlação.
Salario=c(19.307,31.769,22.769,31.307,27.769,30.923,26.508,22.230,28.538,32.307,28.230,19.076,25.384,25.692,42.230,40.923,36.000,47.076,31.461,29.923,47.461,41.153,23.615,40.923,45.076,29.076,44.846)
Experiencia=c(0,17,8,15,9,15,8,5,13,20,11,1,6,7,23,20,18,27,11,10,29,23,4,22,25,9,25)
cbind(Salario,Experiencia)## Salario Experiencia
## [1,] 19.307 0
## [2,] 31.769 17
## [3,] 22.769 8
## [4,] 31.307 15
## [5,] 27.769 9
## [6,] 30.923 15
## [7,] 26.508 8
## [8,] 22.230 5
## [9,] 28.538 13
## [10,] 32.307 20
## [11,] 28.230 11
## [12,] 19.076 1
## [13,] 25.384 6
## [14,] 25.692 7
## [15,] 42.230 23
## [16,] 40.923 20
## [17,] 36.000 18
## [18,] 47.076 27
## [19,] 31.461 11
## [20,] 29.923 10
## [21,] 47.461 29
## [22,] 41.153 23
## [23,] 23.615 4
## [24,] 40.923 22
## [25,] 45.076 25
## [26,] 29.076 9
## [27,] 44.846 25
Uma maneira fácil de obter algumas estatísticas descritivas das variáveis em estudo é através do comando summary(), que retorna as estatísticas mínimo, quartis, média e máximo. Para medir a variabilidade, utilize as funções var() e sd() para obter a variância e o desvio padrão.
summary(Salario)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 19.08 26.10 30.92 32.28 40.92 47.46var(Salario)## [1] 73.68296sd(Salario)## [1] 8.583878Diagrama de Dispersão
plot(Experiencia,Salario, main="Diagrama de Dispersão")
Regressão Linear Simples
Ajuste do Modelo Linear
Sejam X e Y, respectivamente, as variáveis Experiência (explicativa) e Salário (resposta). Propõe-se um modelo de regressão linear de primeira ordem, dado pela equação: \(Y = β_0 + β_1X + Ꞓ\), onde \(β_0\) e \(β_1\) são parâmetros desconhecidos e Ꞓ é o erro aleatório. Para ajustar um modelo de regressão linear no R utiliza-se a função lm:
ajuste=lm(Salario ~ Experiencia)
ajuste##
## Call:
## lm(formula = Salario ~ Experiencia)
##
## Coefficients:
## (Intercept) Experiencia
## 18.061 1.008
Note que função \(lm()\) é chamada com o formato \(lm(y ~ x)\), ou seja, a variável resposta é y e a preditora é x, sempre nessa ordem. O R retorna o valor dos coeficientes de \(β_0\) e \(β_1\) estimados via Método de Mínimos Quadrados. Logo, a equação da reta ajustada é dada por \(Y = 1,83+0,0998Xi\). Com a função summary, diversas medidas descritivas úteis para a análise do ajuste podem ser obtidas:
options(scipen = 999)
summary(ajuste)##
## Call:
## lm(formula = Salario ~ Experiencia)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.9077 -0.8923 0.5774 1.5582 2.7083
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 18.06074 0.81601 22.13 <0.0000000000000002 ***
## Experiencia 1.00770 0.05013 20.10 <0.0000000000000002 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.113 on 25 degrees of freedom
## Multiple R-squared: 0.9417, Adjusted R-squared: 0.9394
## F-statistic: 404 on 1 and 25 DF, p-value: < 0.00000000000000022
Da execução desse comando, pode-se obter, por exemplo, os erros-padrão (Std. Error) das estimativas dos coeficientes de regressão: \(EP(β_0) = 0,0761\) e \(EP(β_1) = 0,0047\). Além disso, obtém-se o valor do Coeficiente de Determinação (Multiple R-Squared), \(R^2 = 0,9478\).Com a função ANOVA, pode-se construir a Tabela da Análise de Variância:
options(scipen = 999)
anova(ajuste)## Analysis of Variance Table
##
## Response: Salario
## Df Sum Sq Mean Sq F value Pr(>F)
## Experiencia 1 1804.12 1804.12 404.01 < 0.00000000000000022 ***
## Residuals 25 111.64 4.47
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Da tabela ANOVA, obtém-se o Quadrado Médio (Mean Sq) Residual, que é uma estimativa para a variância dos erros \((σ^2)\), ou seja, \(s^2 = 0,04\). Para esboçar a reta ajustada no diagrama de dispersão, utilize a função abline:
windows()
plot(Experiencia,Salario)
abline(lm(Salario~Experiencia))
Intervalos de Confiança para βo e β1
Para construir os Intervalos de Confiança (95%) para os coeficientes da regressão, utiliza-se o seguinte comando:
confint(ajuste)## 2.5 % 97.5 %
## (Intercept) 16.3801241 19.741348
## Experiencia 0.9044423 1.110949Análise dos Resíduos
Para avaliar as suposições de que os erros possuem variância constante e são não correlacionados entre si, construa os gráficos de “Resíduos versus Valores Ajustados da Variável Resposta” e “Resíduos versus Valores da Variável Explicativa”:
par(mfrow=c(1, 2))
plot(fitted(ajuste),residuals(ajuste),xlab="Valores Ajustados",ylab="Resíduos")
abline(h=0)
plot(Experiencia,residuals(ajuste),xlab="Experiência",ylab="Resíduos")
abline(h=0)
Para exibir os Valores Ajustados e os Resíduos do ajuste, digite os comandos:
ajuste$residuals ## 1 2 3 4 5
## 1.246263977 -3.422565854 -3.353303002 -1.869174109 0.639001126
## 6 7 8 9 10
## -2.253174109 0.385696998 -0.869215385 -2.622782364 -5.907653471
## 11 12 13 14 15
## -0.915390619 0.007568105 1.277088743 0.577392871 0.992258912
## 16 17 18 19 20
## 2.708346529 -0.199261726 1.807475422 2.315609381 1.785305253
## 21 22 23 24 25
## 0.177083677 -0.084741088 1.523480488 0.692954784 1.822867167
## 26 27
## 1.946001126 1.592867167ajuste$fitted.values## 1 2 3 4 5 6 7 8
## 18.06074 35.19157 26.12230 33.17617 27.13000 33.17617 26.12230 23.09922
## 9 10 11 12 13 14 15 16
## 31.16078 38.21465 29.14539 19.06843 24.10691 25.11461 41.23774 38.21465
## 17 18 19 20 21 22 23 24
## 36.19926 45.26852 29.14539 28.13769 47.28392 41.23774 22.09152 40.23005
## 25 26 27
## 43.25313 27.13000 43.25313Para avaliar a suposição de normalidade dos erros, deve-se construir o gráfico da “Probabilidade Normal dos Resíduos”:
qqnorm(residuals(ajuste), ylab="Resíduos",xlab="Quantis teóricos",main="")
qqline(residuals(ajuste))
Pela figura acima (página seguinte), observa-se a violação da suposição de que os erros aleatórios têm distribuição Normal. Considere, também o Teste de Normalidade de Shapiro Wilk:
shapiro.test(residuals(ajuste))##
## Shapiro-Wilk normality test
##
## data: residuals(ajuste)
## W = 0.90144, p-value = 0.01444Portanto, como o Valor P do teste é pequeno, rejeita-se a hipótese de normalidade dos resíduos e, por consequência, conclui-se que os erros não são normalmente distribuídos.