Aula de Modelos Lineares I
Patrícia Lusié Velozo da Costa
2025-07-02
Modelo de Regressão Linear Múltipla - Medidas Corretivas para heterocedasticidade
1 Simulando um conjunto de dados heterocedástico
Simule um conjunto de dados considerando \[Y_i \stackrel{indep}{\sim} N(\beta_0 + \beta_1 X_{i1} + \beta_2 X_{i2}, \sigma^2 X_{i2}^2), i=1, \ldots n,\] sendo n=1000 (tamanho amostral), \(\beta_0=2\), \(\beta_1=5\), \(\beta_2=3\), \(\sigma^2=4\), \(X_{i1} \stackrel{iid}{\sim} N(0,1)\) e \(X_{i2} \stackrel{iid}{\sim} t_{(5)}\).
# Configurações iniciais
n = 1000 # Número de observações
# Simulando as covariáveis
x1 = rnorm(n, mean = 0, sd = sqrt(1)) # Covariável 1
x2 = rt(n, 5) # Covariável 2
x = cbind(1, x1=x1, x2=x2)
#definindo total de covariaveis
p = ncol(x)-1
# Definindo os coeficientes reais
beta_0 = 2 # Intercepto
beta_1 = 5 # Coeficiente para x1
beta_2 = 3 # Coeficiente para x2
beta = c(beta_0, beta_1, beta_2)
# Gerando a variável resposta com erro aleatório
sigma2 = 4
e = rnorm(n, mean = 0, sd = sqrt(sigma2*x2^2)) # Erro aleatório
y = x %*% beta + e
# Criando um data frame com os dados
dados = data.frame(y, x1=x1, x2=x2)2 Analisando os dados
Analise o conjunto de dados simulado.
#analisando os dados gerados
par(mfrow=c(1,2))
hist(y, main="", freq=F, bty="n", ylab="densidade")
hist(e, main="", freq=F, bty="n", ylab="densidade")
par(mfrow=c(1,3))
plot(x1, y, bty="n")
plot(x2, y, bty="n")
plot(x1, x2, bty="n")
round(cor(dados),3)## y x1 x2
## y 1.000 0.739 0.572
## x1 0.739 1.000 -0.013
## x2 0.572 -0.013 1.000shapiro.test(e)##
## Shapiro-Wilk normality test
##
## data: e
## W = 0.90087, p-value < 2.2e-16ks.test(e,"pnorm",0,sqrt(sigma2))##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: e
## D = 0.10893, p-value = 9.856e-11
## alternative hypothesis: two-sided3 Ajustando o MRLM usando inferência clássica
Suponha que o modelo usado para a geração dos dados seja desconhecido. Ajuste o seguinte modelo de regressão linear múltipla \[Y_i \stackrel{ind}{\sim} N(\gamma_0+\gamma_1 X_{i1}+\gamma_2 X_{i2}, \delta^2),\] sendo \(\gamma_0, \gamma_1, \gamma_2, \delta^2\) desconhecidos. Em seguida, faça uma análise dos resíduos.
1. Será que os coeficientes ajustados \((\hat{\gamma}_0, \hat{\gamma}_1, \hat{\gamma}_2)\) e a variância estimada \(\hat{\delta}^2\) são bons estimadores para os efeitos \((\beta_0, \beta_1, \beta_2)\) e a variância \(\sigma^2\) do modelo que originou os dados?
2. Sejam \(\hat{\gamma}_0, \hat{\gamma}_1, \hat{\gamma}_2\) estimadores de mínimos quadrados para \(\gamma_0, \gamma_1, \gamma_2\), respectivamente. Você espera encontrar homocedasticidade nos resíduos \(\epsilon_i = Y_i - \hat{\gamma}_0 - \hat{\gamma}_1 X_{i1} - \hat{\gamma}_2 X_{i2}\) deste ajuste?
#ajustando o modelo
ajuste = lm(y ~ x1 + x2, data = dados)
summary(ajuste)##
## Call:
## lm(formula = y ~ x1 + x2, data = dados)
##
## Residuals:
## Min 1Q Median 3Q Max
## -13.9460 -0.7834 -0.0402 0.7640 12.1427
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.00544 0.07141 28.08 <2e-16 ***
## x1 5.05345 0.07284 69.38 <2e-16 ***
## x2 2.94581 0.05440 54.15 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.258 on 997 degrees of freedom
## Multiple R-squared: 0.8847, Adjusted R-squared: 0.8844
## F-statistic: 3823 on 2 and 997 DF, p-value: < 2.2e-16#beta do modelo que originou os dados versus gamachapeu e seu IC
alfa = 0.05
cbind(beta, ajuste$coefficients, confint(ajuste, level=1-alfa))## beta 2.5 % 97.5 %
## (Intercept) 2 2.005440 1.865307 2.145573
## x1 5 5.053445 4.910515 5.196376
## x2 3 2.945811 2.839061 3.052562#sigma2 do modelo que originou os dados versus delta2chapeu e seu IC
c(sigma2, mean(sigma2*x2^2), summary(ajuste)$sigma^2,
(n-p-1)*summary(ajuste)$sigma^2/qchisq(1-alfa/2, n-p-1),
(n-p-1)*summary(ajuste)$sigma^2/qchisq(alfa/2, n-p-1))## [1] 4.000000 6.894087 5.098982 4.679389 5.577954#Analisando os resíduos
hist(ajuste$residuals, main="", freq=F, ylab="densidade", xlab="resíduos")
curve(dnorm(x, 0, sqrt(summary(ajuste)$sigma^2)), add=T, col="red", lwd=2)
qqnorm(ajuste$residuals, xlab="quantis teóricos", ylab="quantis amostrais", bty="n")
qqline(ajuste$residuals, col="red", lwd=2)
shapiro.test(ajuste$residuals)##
## Shapiro-Wilk normality test
##
## data: ajuste$residuals
## W = 0.90027, p-value < 2.2e-16ks.test(ajuste$residuals,"pnorm",0, sqrt(summary(ajuste)$sigma^2))##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: ajuste$residuals
## D = 0.12359, p-value = 1.081e-13
## alternative hypothesis: two-sidedpar(mfrow=c(1,4))
plot(x1, ajuste$residuals)
abline(h=0, lwd=2, col="red", lty=3)
plot(x2, ajuste$residuals)
abline(h=0, lwd=2, col="red", lty=3)
plot(y, ajuste$residuals)
abline(h=0, lwd=2, col="red", lty=3)
plot(ajuste$residuals)
abline(h=0, lwd=2, col="red", lty=3)
# Testando Homocedasticidade (Teste de Breusch-Pagan)
library(lmtest)
bptest(ajuste)##
## studentized Breusch-Pagan test
##
## data: ajuste
## BP = 19.106, df = 2, p-value = 7.1e-05#Interpretação: Um p-valor pequeno indica a presença de heterocedasticidade, o que significa que a variância dos resíduos não é constante.
#Avaliando a Independência dos Resíduos (Durbin-Watson)
dwtest(ajuste)##
## Durbin-Watson test
##
## data: ajuste
## DW = 1.9391, p-value = 0.1673
## alternative hypothesis: true autocorrelation is greater than 0#Interpretação: Um valor de estatística Durbin-Watson próximo de 2 sugere que não há autocorrelação dos resíduos.4 Ajustando o MRLM heterocedástico
Suponha agora que a variância seja conhecida e que deseja-se então ajustar o seguinte modelo de regressão linear múltipla \[Y_i = \beta_0+\beta_1 X_{i1}+\beta_2 X_{i2} + e_i, \quad e_i \stackrel{ind}{\sim}N(0, \sigma^2 X_{i2}^2),\] sendo \(\beta_0, \beta_1, \beta_2\) desconhecidos. Ajuste isso e faça uma análise de resíduos.
OBS:
Supondo \(V_i = \sigma^2 X_{i2}^2\) conhecido, tem-se que \[R_i = \frac{e_i}{\sqrt{V_i}} \stackrel{iid}{\sim}N(0, 1).\]
Note que os resíduos do modelo ajustado serão iguais a \[\hat{e}_i = Y_i - \hat{\beta}_0-\hat{\beta}_1 X_{i1} -\hat{\beta}_2 X_{i2}.\] Considere \[\hat{R}_i = \frac{\hat{e}_i}{\sqrt{V_i}} = \frac{Y_i - \hat{\beta}_0-\hat{\beta}_1 X_{i1}-\hat{\beta}_2 X_{i2}}{\sqrt{V_i}}.\]
Como \(E[\hat{\beta}_j] = \beta_j\), \(\forall j=0, \ldots, p\), tem-se que \[E[\hat{R}_i] = \frac{Y_i - \beta_0-\beta_1 X_{i1} -\beta_2 X_{i2}}{\sqrt{V_i}} = \frac{e_i}{\sqrt{V_i}} = R_i.\] Logo, para verificar a normalidade, a homocedasticidade e a independência, avalie os resíduos modificados \(\hat{R}_i\).
#ajustando o modelo
omegai = 1/(sigma2*x2^2)
ajusteHV = lm(y ~ x1 + x2, data = dados, weights = omegai)
#analisando o ajuste
summary(ajusteHV)##
## Call:
## lm(formula = y ~ x1 + x2, data = dados, weights = omegai)
##
## Weighted Residuals:
## Min 1Q Median 3Q Max
## -3.2127 -0.7061 -0.0446 0.6490 3.2080
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.999979 0.001448 1381.37 <2e-16 ***
## x1 5.000215 0.002085 2397.66 <2e-16 ***
## x2 2.994906 0.061432 48.75 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9711 on 997 degrees of freedom
## Multiple R-squared: 0.9998, Adjusted R-squared: 0.9998
## F-statistic: 2.875e+06 on 2 and 997 DF, p-value: < 2.2e-16#Note que o modelo ajustado foi o modelo que usamos para gerar os dados.
#beta versus betachapeu e seu IC
alfa = 0.05
cbind(beta, ajusteHV$coefficients, confint(ajusteHV, level=1-alfa))## beta 2.5 % 97.5 %
## (Intercept) 2 1.999979 1.997138 2.002820
## x1 5 5.000215 4.996123 5.004308
## x2 3 2.994906 2.874354 3.115457#esperamos que a variancia estimada seja 1
#sigma2chapeu e seu IC
c(1, summary(ajusteHV)$sigma^2,
(n-p-1)*summary(ajusteHV)$sigma^2/qchisq(1-alfa/2, n-p-1),
(n-p-1)*summary(ajusteHV)$sigma^2/qchisq(alfa/2, n-p-1))## [1] 1.0000000 0.9431063 0.8654985 1.0316969#analisando os residuos homocedasticos
residuos_Richapeu = ajusteHV$residuals*sqrt(omegai)
# Analisando os resíduos
hist(residuos_Richapeu, main="", freq=F, ylab="densidade", xlab="resíduos")
curve(dnorm(x, 0, sqrt(summary(ajusteHV)$sigma^2)), add=T, col="red", lwd=2)
qqnorm(residuos_Richapeu, xlab="quantis teóricos", ylab="quantis amostrais", bty="n")
qqline(residuos_Richapeu, col="red", lwd=2)
shapiro.test(residuos_Richapeu)##
## Shapiro-Wilk normality test
##
## data: residuos_Richapeu
## W = 0.99855, p-value = 0.5858ks.test(residuos_Richapeu,"pnorm",0, sqrt(summary(ajusteHV)$sigma^2))##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: residuos_Richapeu
## D = 0.026754, p-value = 0.4713
## alternative hypothesis: two-sidedpar(mfrow=c(1,4))
plot(x1, residuos_Richapeu)
abline(h=0, lwd=2, col="red", lty=3)
plot(x2, residuos_Richapeu)
abline(h=0, lwd=2, col="red", lty=3)
plot(y, residuos_Richapeu)
abline(h=0, lwd=2, col="red", lty=3)
plot(residuos_Richapeu)
abline(h=0, lwd=2, col="red", lty=3)
5 Ajustando o MRLM heterocedástico - forma 2
Ajustar o modelo de regressão linear múltipla \[Y_i = \beta_0+\beta_1 X_{i1} +\beta_2 X_{i2} + e_i, \quad e_i \stackrel{ind}{\sim}N(0, V_i = \sigma^2 X_{i2}^2)\] equivale a ajustar o modelo \[\underbrace{\frac{Y_i}{\sqrt{V_i}}}_{Y_i^{*}} = \beta_0 \underbrace{\frac{1}{\sqrt{V_i}}}_{X_{i0}^*} + \beta_1 \underbrace{ \frac{X_{i1}}{\sqrt{V_i}} }_{X_{i1}^*} +\beta_2 \underbrace{\frac{X_{i2}}{\sqrt{V_i}}}_{X_{i2}^*} + \underbrace{\frac{e_i}{\sqrt{V_i}}}_{e_i^*}\] \[ Y_i^{*} = \beta_0 X_{i0}^* + \beta_1 X_{i1}^* +\beta_2 X_{i2}^* + e_i^*, \quad e_i^* \stackrel{ind}{\sim}N(0, 1).\]
#ajustando o modelo
y_mod = sqrt(omegai)*y
x_mod = sqrt(omegai)*x
dados_mod = data.frame(y = y_mod, intercepto = sqrt(omegai), x1 = x_mod[,2], x2 = x_mod[,3])
ajusteHV_mod = lm(y ~ intercepto + x1 + x2 - 1, data = dados_mod)
#analisando o ajuste
summary(ajusteHV_mod)##
## Call:
## lm(formula = y ~ intercepto + x1 + x2 - 1, data = dados_mod)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.2127 -0.7061 -0.0446 0.6490 3.2080
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## intercepto 1.999979 0.001448 1381.37 <2e-16 ***
## x1 5.000215 0.002085 2397.66 <2e-16 ***
## x2 2.994906 0.061432 48.75 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9711 on 997 degrees of freedom
## Multiple R-squared: 1, Adjusted R-squared: 1
## F-statistic: 8.7e+06 on 3 and 997 DF, p-value: < 2.2e-16#Note que o modelo ajustado foi o modelo que usamos para gerar os dados.
#beta versus betachapeu e seu IC
alfa = 0.05
cbind(beta, ajusteHV$coefficients, confint(ajusteHV, level=1-alfa), ajusteHV_mod$coefficients, confint(ajusteHV_mod, level=1-alfa))## beta 2.5 % 97.5 % 2.5 % 97.5 %
## (Intercept) 2 1.999979 1.997138 2.002820 1.999979 1.997138 2.002820
## x1 5 5.000215 4.996123 5.004308 5.000215 4.996123 5.004308
## x2 3 2.994906 2.874354 3.115457 2.994906 2.874354 3.115457#esperamos que a variancia estimada seja 1
#sigma2chapeu e seu IC
c(1, summary(ajusteHV)$sigma^2,
(n-p-1)*summary(ajusteHV)$sigma^2/qchisq(1-alfa/2, n-p-1),
(n-p-1)*summary(ajusteHV)$sigma^2/qchisq(alfa/2, n-p-1),
summary(ajusteHV_mod)$sigma^2,
(n-p-1)*summary(ajusteHV_mod)$sigma^2/qchisq(1-alfa/2, n-p-1),
(n-p-1)*summary(ajusteHV_mod)$sigma^2/qchisq(alfa/2, n-p-1))## [1] 1.0000000 0.9431063 0.8654985 1.0316969 0.9431063 0.8654985 1.0316969#analisando os residuos homocedasticos
residuos_eichapeu = ajusteHV_mod$residuals
# Analisando os resíduos
hist(residuos_eichapeu, main="", freq=F, ylab="densidade", xlab="resíduos")
curve(dnorm(x, 0, sqrt(summary(ajusteHV_mod)$sigma^2)), add=T, col="red", lwd=2)
qqnorm(residuos_eichapeu, xlab="quantis teóricos", ylab="quantis amostrais", bty="n")
qqline(residuos_eichapeu, col="red", lwd=2)
shapiro.test(residuos_eichapeu)##
## Shapiro-Wilk normality test
##
## data: residuos_eichapeu
## W = 0.99855, p-value = 0.5858ks.test(residuos_eichapeu,"pnorm",0, sqrt(summary(ajusteHV_mod)$sigma^2))##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: residuos_eichapeu
## D = 0.026754, p-value = 0.4713
## alternative hypothesis: two-sidedpar(mfrow=c(1,4))
plot(x1, residuos_eichapeu)
abline(h=0, lwd=2, col="red", lty=3)
plot(x2, residuos_eichapeu)
abline(h=0, lwd=2, col="red", lty=3)
plot(y, residuos_eichapeu)
abline(h=0, lwd=2, col="red", lty=3)
plot(residuos_eichapeu)
abline(h=0, lwd=2, col="red", lty=3)
# Testando Homocedasticidade (Teste de Breusch-Pagan)
bptest(ajusteHV_mod)##
## studentized Breusch-Pagan test
##
## data: ajusteHV_mod
## BP = 3.1873, df = 2, p-value = 0.2032#Interpretação: Um p-valor pequeno indica a presença de heterocedasticidade, o que significa que a variância dos resíduos não é constante.
#Avaliando a Independência dos Resíduos (Durbin-Watson)
dwtest(ajusteHV_mod)##
## Durbin-Watson test
##
## data: ajusteHV_mod
## DW = 1.8591, p-value = 0.01388
## alternative hypothesis: true autocorrelation is greater than 0#Interpretação: Um valor de estatística Durbin-Watson próximo de 2 sugere que não há autocorrelação dos resíduos.6 Ajustando o MRLM heterocedástico supondo sigma2 desconhecido
Suponha agora que deseja-se ajustar o seguinte modelo de regressão linear múltipla \[Y_i = \beta_0+\beta_1 X_{i1}+\beta_2 X_{i2} + e_i, \quad e_i \stackrel{ind}{\sim}N(0, \sigma^2 X_{i2}^2),\] sendo \(\beta_0, \beta_1, \beta_2, \sigma^2\) desconhecidos. Ajuste isso e faça uma análise de resíduos.
OBS: Tem-se que \[R_i^* = \frac{e_i}{\sqrt{X_{i2}^2}} \stackrel{iid}{\sim}N(0, \sigma^2).\] Logo, os efeitos \(R_i^*\) são homocedásticos, independentes e normalmente distribuídos.
#ajustando o modelo
omegai = 1/(x2^2)
ajusteHVs = lm(y ~ x1 + x2, data = dados, weights = omegai)
#analisando o ajuste
summary(ajusteHVs)##
## Call:
## lm(formula = y ~ x1 + x2, data = dados, weights = omegai)
##
## Weighted Residuals:
## Min 1Q Median 3Q Max
## -6.4254 -1.4121 -0.0893 1.2981 6.4160
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.999979 0.001448 1381.37 <2e-16 ***
## x1 5.000215 0.002085 2397.66 <2e-16 ***
## x2 2.994906 0.061432 48.75 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.942 on 997 degrees of freedom
## Multiple R-squared: 0.9998, Adjusted R-squared: 0.9998
## F-statistic: 2.875e+06 on 2 and 997 DF, p-value: < 2.2e-16#beta versus betachapeu e seu IC
alfa = 0.05
cbind(beta, ajusteHVs$coefficients, confint(ajusteHVs, level=1-alfa))## beta 2.5 % 97.5 %
## (Intercept) 2 1.999979 1.997138 2.002820
## x1 5 5.000215 4.996123 5.004308
## x2 3 2.994906 2.874354 3.115457#esperamos que a variancia estimada seja sigma2
#sigma2chapeu e seu IC
c(sigma2, summary(ajusteHVs)$sigma^2,
(n-p-1)*summary(ajusteHVs)$sigma^2/qchisq(1-alfa/2, n-p-1),
(n-p-1)*summary(ajusteHVs)$sigma^2/qchisq(alfa/2, n-p-1))## [1] 4.000000 3.772425 3.461994 4.126788#analisando os residuos homocedasticos
residuos_Richapeu_s = ajusteHVs$residuals*sqrt(omegai)
# Analisando os resíduos
hist(residuos_Richapeu_s, main="", freq=F, ylab="densidade", xlab="resíduos")
curve(dnorm(x, 0, sqrt(summary(ajusteHVs)$sigma^2)), add=T, col="red", lwd=2)
qqnorm(residuos_Richapeu_s, xlab="quantis teóricos", ylab="quantis amostrais", bty="n")
qqline(residuos_Richapeu_s, col="red", lwd=2)
shapiro.test(residuos_Richapeu_s)##
## Shapiro-Wilk normality test
##
## data: residuos_Richapeu_s
## W = 0.99855, p-value = 0.5858ks.test(residuos_Richapeu_s,"pnorm",0, sqrt(summary(ajusteHVs)$sigma^2))##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: residuos_Richapeu_s
## D = 0.026754, p-value = 0.4713
## alternative hypothesis: two-sidedpar(mfrow=c(1,4))
plot(x1, residuos_Richapeu_s)
abline(h=0, lwd=2, col="red", lty=3)
plot(x2, residuos_Richapeu_s)
abline(h=0, lwd=2, col="red", lty=3)
plot(y, residuos_Richapeu_s)
abline(h=0, lwd=2, col="red", lty=3)
plot(residuos_Richapeu_s)
abline(h=0, lwd=2, col="red", lty=3)
7 Ajustando o MRLM com a correção de heterocedasticidade
Do modelo \[Y_i \stackrel{ind}{\sim} N(\gamma_0+\gamma_1 X_{i1}+\gamma_2 X_{i2}, \delta^2)\] obtenha os resíduos \(\hat{\zeta}_i = Y_i - \hat{Y}_i^*\) sendo \(\hat{Y}_i^* = \hat{\gamma}_0+\hat{\gamma}_1 X_{i1}+\hat{\gamma}_2 X_{i2}\).
Em seguida, ajuste o modelo \[|\hat{\zeta}_i| \stackrel{ind}{\sim} N(\alpha_0+\alpha_1 \hat{Y}_i^*, \eta^2).\] Defina \(\hat{s}_i = \hat{\alpha}_0 + \hat{\alpha}_1 \hat{Y}_i^*\).
Defina \(\omega_i = \frac{1}{\hat{s}_i^2}\) e ajuste o modelo \[Y_i = \beta_0+\beta_1 X_{i1}+\beta_2 X_{i2} + e_i, \quad e_i \stackrel{ind}{\sim}N\left(0, \frac{1}{\omega_i}\right).\]
#ajustando o modelo para os residuos do modelo ajustado considerando homocedasticidade
ajusteH = lm(abs(ajuste$residuals) ~ ajuste$fitted.values, data = dados)
si = ajusteH$fitted.values
omega_i = 1/(si^2)
ajusteH_2 = lm(y ~ x1+x2, data = dados, weights = omega_i)
summary(ajusteH_2)##
## Call:
## lm(formula = y ~ x1 + x2, data = dados, weights = omega_i)
##
## Weighted Residuals:
## Min 1Q Median 3Q Max
## -9.8145 -0.5510 -0.0283 0.5377 8.5489
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.00544 0.07141 28.08 <2e-16 ***
## x1 5.05336 0.07284 69.37 <2e-16 ***
## x2 2.94583 0.05441 54.14 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.589 on 997 degrees of freedom
## Multiple R-squared: 0.8846, Adjusted R-squared: 0.8844
## F-statistic: 3822 on 2 and 997 DF, p-value: < 2.2e-16#analisando os residuos homocedasticos
residuosH_Richapeu = ajusteH_2$residuals * sqrt(omegai)
hist(residuosH_Richapeu, main="", freq=F, ylab="densidade", xlab="resíduos")
curve(dnorm(x, 0, sqrt(summary(ajusteH_2)$sigma^2)), add=T, col="red", lwd=2)
qqnorm(residuosH_Richapeu, xlab="quantis teóricos", ylab="quantis amostrais", bty="n")
qqline(residuosH_Richapeu, col="red", lwd=2)
shapiro.test(residuosH_Richapeu)##
## Shapiro-Wilk normality test
##
## data: residuosH_Richapeu
## W = 0.55548, p-value < 2.2e-16ks.test(residuosH_Richapeu,"pnorm",0, sqrt(summary(ajusteH_2)$sigma^2))##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: residuosH_Richapeu
## D = 0.079928, p-value = 5.651e-06
## alternative hypothesis: two-sidedpar(mfrow=c(1,4))
plot(x1, residuosH_Richapeu)
abline(h=0, lwd=2, col="red", lty=3)
plot(x2, residuosH_Richapeu)
abline(h=0, lwd=2, col="red", lty=3)
plot(y, residuosH_Richapeu)
abline(h=0, lwd=2, col="red", lty=3)
plot(residuosH_Richapeu)
abline(h=0, lwd=2, col="red", lty=3)