# Aula prática de SME 260 - Análise de Regressão 04/04/2011
# Profa. Cibele Russo

# Faça download do arquivo AulaPrática.zip na WIKI 

# Extraia os arquivos para o diretório Z: ou outro de sua preferência



# ------------- Item 1 --------------
 

dados<-read.table("D:\\AulaPrática\\censo.txt", header=TRUE)

attach(dados)

Y <- renda

X <- escolaridade

n <- length(Y)

# Ajuste do modelo

 fit.model <- lm(Y~X)

 summary(fit.model)


# Estimativa de sigma

 sigma.chapéu <- sqrt(sum(fit.model$residuals^2)/(n-2))
 
 plot(X,Y, pch=16) 
 
 abline(fit.model) # Reta ajustada pelo R com a função lm: linear model
 

# Tabela ANOVA

 anova(fit.model) 

 
# Falta de ajuste e erro puro

 # install.packages("alr3") ## Executar se necessário para instalar a biblioteca 
 
 library(alr3)
 
 pureErrorAnova(fit.model) 
 

# Análise de resíduos

 y.ajust <- fitted(fit.model)

 resíduos <- residuals(fit.model)

 plot(y.ajust, resíduos)

 plot(y.ajust, resíduos, pch=16, xlab=expression(hat(Y)), ylab=expression(e[i]), main="Gráfico de resíduos")

 rp <- resíduos/sigma.chapéu # resíduos padronizados ordenados

 faixay <- range(c(-3,3),rp)

 plot(y.ajust, rp, pch=16, ylim=faixay, xlab=expression(hat(Y)), ylab=expression(d[i]), main="Gráfico de resíduos padronizados")

  abline(qnorm(0.975),0, col=2, lty=2)

  abline(qnorm(0.025),0, col=2, lty=2)

  #  identify(y.ajust,rp)


# Gráfico de probabilidade normal 

 rp<-sort(rp)

 q <- qnorm( (1:n - 0.5)/n ) # quantis teóricos da normal padrão

 plot(q, rp, pch=16, main="Gráfico de probabilidade normal", xlab="Quantis teóricos", ylab="Resíduos padronizados")

 abline(0,1)

 # qqnorm(rp, pch=16)   # faz direto o gráfico QQ 

# Envelope com resíduos estudentizados

 source("D:\\AulaPrática\\envel_norm.txt") 
 



# ------------- Item 2 --------------

 rm(list=ls(all=TRUE)) #remove todos os objetos existentes 

 dados <- read.table("D:\\AulaPrática\\tratores.txt", header=TRUE)

 attach(dados)
 
 Y <- custo

 X <- idade

 n <- length(Y)


# ------------- Item 3 --------------

 rm(list=ls(all=TRUE)) #remove todos os objetos existentes 

 dados <- read.table("D:\\AulaPrática\\foguetes.txt", header=TRUE)

 attach(dados)
 
 Y <- resistência

 X <- idade

 n <- length(Y)



# ------------- Item 4 --------------

rm(list=ls(all=TRUE)) #remove todos os objetos existentes 

 Acuidade <- read.table("D:\\AulaPrática\\Acuidade.txt", header=TRUE)

 # Alternativa para leitura de dados separados por vírgulas
 # Acuidade <- read.csv("D:\\AulaPrática\\Acuidade.csv", header=TRUE, sep=";")

 attach(Acuidade)

 # attach(Acuidade[which(sexo=="M"),]) # Para considerar somente dados de um sexo
 # attach(Acuidade[which(sexo=="F"),])

 Y <- tempo

 X <- idade

 n <- length(Y)