#***********************************************************************
# Laboratorio 14.  Minimos cuadrados ponderados
#  Edgar Acuna, Marzo 2003
#  
#*********************************************************************
#Leyendo el archivo de datos de millaje directamente de la internet
millaje<-read.table(file="http://math.uprm.edu/~edgar/millaje.txt",header=T)
#millaje<-read.table(file="c:/millaje.txt",header=T)
#
#Extrayendo la primera y ultima observacion y consideramdo solo
#las columnas de variable de respuesta(1) y de wt(3).
#
millaje1<-millaje[-c(1,82),c(1,3)]
win.graph()
plot(millaje1[,2],millaje1[,1])
#Agrupando los valores de x
tabla<-table(millaje1[,2])
valsx<-dimnames(tabla)
valsx<-as.double(unlist(valsx))
nvalsx<-length(valsx)
#Calculando la varianza de cada grupo
varsg<-rep(0,nvalsx)
for(i in 1:nvalsx)
{grupox<-millaje1[millaje1[,2]==valsx[i],]
varsg[i]<-var(grupox[,1])
}
varsg
#Asignando los pesos a cada observacion (metodo I)
win.graph()
plot(varsg,valsx)
ndatos<-dim(millaje1)[1]
pesos<-rep(0,ndatos)
for(i in 1:ndatos)
{for(j in 1:nvalsx)
{if(millaje1[i,2]==valsx[j])
{pesos[i]<-1/varsg[j]}
}
}
pesos
#resultados de regresion sin ponderar
lw1<-lm(mpg~.,data=millaje1,weights=pesos)
summary(lw1)
#Asignando los pesos a cada observacion (metodo II)
valsx2<-valsx^2
s2<-lm(varsg~valsx+valsx2)
s2
valsx<-millaje1[,2]
valsx2<-millaje1[,2]^2
nuevo<-data.frame(valsx,valsx2)
pesos1<-predict.lm(s2,nuevo)
pesos1<-1/pesos1
pesos1
#resultados de regresion ponderada
lw2<-lm(mpg~.,data=millaje1,weights=pesos1)
summary(lw2)
#************************************************
#Aplicando regresion ponderada con la variable sp
#************************************************
nclases<-nclass.scott(millaje[,2])
clases<-cut(millaje[,2],nclases)
millaje2<-cbind(millaje[,1],clases)
#Calculando la varianza de cada grupo
varsg2<-rep(0,nclases)
for(i in 1:nclases)
{grupox<-millaje2[millaje2[,2]==i,]
varsg2[i]<-var(grupox[,1])
}
varsg2
#Asignando los pesos a cada observacion 
win.graph()
plot(varsg2,1:nclases)
nd<-dim(millaje)[1]
pesos2<-rep(0,nd)
for(i in 1:nd)
{for(j in 1:nclases)
{if(millaje2[i,2]==j)
{pesos2[i]<-1/varsg2[j]}
}
}
pesos2
#Resultados usando regresion ponderada
lw3<-lm(mpg~sp,data=millaje[,c(1,2)],weights=pesos2)
summary(lw3)
#Resultados sin usar regresion ponderada
lw4<-lm(mpg~sp,data=millaje)
summary(lw4)