Não foi possível enviar o arquivo. Será algum problema com as permissões?
Diferenças

Diferenças

Aqui você vê as diferenças entre duas revisões dessa página.

Link para esta página de comparações

cursos:rbelem:dia3a [2010/05/27 11:59]
paulojus criada
cursos:rbelem:dia3a [2010/05/29 11:20] (atual)
paulojus
Linha 1: Linha 1:
 +====== Dia 4 - a ======
 +
 +Vamos revisar, comentar e alterar um código enviado na lista brasileira do do R 
 +para ilustrar alguns aspectos do uso adequado da linguagem.
 +
 +Aproveitamos o exemplo para discutir a especificação de modelos através do uso de fórmulas no R.
 +
 <code R> <code R>
 +## definindo número de dados  ​
 idados<- 50 idados<- 50
-adados<​-20000+##adados<​-20000
    
-#Equação verdadeira de relacionamento (usada para simular Fase I e gerar Modelo de Regressão) +# Equação verdadeira de relacionamento ​ 
-coef0 <- c(10, 2, -3, 2.3)+## deseja simular de y = beta_0 + beta_1 X1 + beta_2 X2 + beta_3 X1.X2  
 +(usada para simular Fase I e gerar Modelo de Regressão) 
 +coef0 <- c(10, 2, -3, 2.3)    ## parâmetros para simular o modelo
    
 #Equação alterada #Equação alterada
-coefal<​-c(10,​ 6, -3, 2.3)+#coefal<​-c(10,​ 6, -3, 2.3)
    
 #Fase I #Fase I
 +## simulando o valor das covariáveis X1 e X2
 x1 <​-runif(idados,​ min=-1, max=1) x1 <​-runif(idados,​ min=-1, max=1)
 x2 <​-runif(idados,​ min=-1, max=1) x2 <​-runif(idados,​ min=-1, max=1)
-  ​+## simulando os erros:  ​
 e <​-rnorm(idados,​ mean=0, sd=4) e <​-rnorm(idados,​ mean=0, sd=4)
  
 +## criando a matriz X 
 x <​-cbind(1,​ x1, x2, x1*x2) x <​-cbind(1,​ x1, x2, x1*x2)
 +## Comentários:​  
 +##     não precisa criar com a coluna x1*x2  
 +##     pode usar a função model.matrix() para montar a matriz a partir 
 +##         da espeficicação do modelo
 X <- model.matrix(~ x1*x2) X <- model.matrix(~ x1*x2)
 X X
  
-model.matrix(~ I(x1*x2))+## Sobre formulas no R: 
 +## Notar que os operadores +, * , ^ etc, tem um sentido diferente em formulas 
 +## I()  "​força"​ a operação algébrica entre as variáveis 
 +head(model.matrix(~ x1*x2)) 
 +head(model.matrix(~ I(x1*x2))) 
 +head(model.matrix(~ x1+x2)) 
 +head(model.matrix(~ I(x1+x2)))
  
  
Linha 26: Linha 47:
 da <​-data.frame (x1, x2, x1*x2) da <​-data.frame (x1, x2, x1*x2)
  
 +## colocando as covariáveis num data-frame
 da <​-data.frame (x1=x1, x1=x2) da <​-data.frame (x1=x1, x1=x2)
 +## Com.:note que poderi já ter criado aqui para evitar objetos desnecessários:​
 +da <​-data.frame(x1=runif(idados,​ min=-1, max=1), ​
 +                 ​x2=runif(idados,​ min=-1, max=1)
 +               )
  
 +## simulando os dados
 y  <- vector(length=idados) y  <- vector(length=idados)
 ii<​-1:​idados ii<​-1:​idados
Linha 34: Linha 61:
 } }
  
-yy <- as.vector(X %*% coef0 + e+## Com: note que não precisa fazer o loop for() pois o R é linguagem vetorial 
-all.equal(y,​ yy) +##       alem disto poderia criar o y direto no data-frame 
 +da$y <- X %*% coef0 + e
 da$y  <-y da$y  <-y
-da$y <- as.vector(X %*% coef0 + e) 
  
-rm(x1, x2, y, yy)+all.equal(yda$y)
  
 +## Com: procure sempre apagar o que não é mais necessário
 +rm(x1, x2, y)
 +
 +## note a diferença (e semelhancas) entre os modelos
 mod   <​-lm(y~ x1 + x2 + x1*x2, data=da) mod   <​-lm(y~ x1 + x2 + x1*x2, data=da)
 +mod
 +lm(y~ x1 + x2 + I(x1*x2), data=da)
 +lm(y~ x1 + x2 + x1:x2, data=da)
 +lm(y~ x1*x2, data=da)
 +lm(y~ (x1 + x2)^2, data=da)
  
-mod   <​-lm(y~ x1 + x2 + I(x1*x2), data=da) 
-mod   <​-lm(y~ x1 + x2 + x1:x2, data=da) 
-mod   <​-lm(y~ x1*x2, data=da) 
  
 sdmod <​-sqrt(sum(residuals(mod)^2)/​df.residual(mod)) sdmod <​-sqrt(sum(residuals(mod)^2)/​df.residual(mod))
 +## .. ou direto de summary()
 sdmod <- summary(mod)$sigma sdmod <- summary(mod)$sigma
   ​   ​
Linha 75: Linha 107:
    
    
-## +## Obtendo as estimativas de beta sem usar lm() usando contas matriciais 
 +## beta = (X'​X)^(-1) X'y 
 +## operação correta porém ​ ineficiente e inadequada:
 solve(t(X) %*% X) %*% t(X) %*% y solve(t(X) %*% X) %*% t(X) %*% y
 +## operação feita de maneira mais adequada:
 solve(crossprod(X),​ crossprod(X,​y)) solve(crossprod(X),​ crossprod(X,​y))
 +</​code>​
  
-</​code>​ 

QR Code
QR Code cursos:rbelem:dia3a (generated for current page)