#***********************************************
#Laboratorio 21: Suavizadores bidimensionales
# regresorgrama, running means, regresion por kernels y lowess
#Edgar Acuna, Mayo 2003
#**********************************************
regresorg<-function (x,y,k) 
{# ***********************************************************
# Funcion que calcula el regresorgrama
# Input 
# Edgar Acuna, Mayo 2003
# ************************************************************
n<-length(x)
x<-sort(x)
y<-y[order(x)]
xpoints<-x[1]
nint<-floor(n/k)
ymeans<-rep(0,k)
for(j in 1:k)
{
ind<-((j-1)*nint+1):(j*nint)
if(j<k)
xpoints<-c(xpoints,x[j*nint])
if(j==k)
{ind<-((j-1)*nint+1):n
xpoints<-c(xpoints,x[n])
}
ymeans[j]<-mean(y[ind])
}
xpoints<-c(xpoints,xpoints[2:k])
xpoints<-sort(xpoints)
ymeans1<-rep(ymeans,each=2)
plot(x,y)
lines(xpoints,ymeans1)
title("Regresorgrama")
cat("\las medias de y en cada subintervalo son:\n")
ymeans
}
#Aplicando regresorgrama al conjunto de datos air
air<-read.table("http://math.uprm.edu/~edgar/air.txt",header=T)
#air<-read.table("c:/air.txt",header=T)
regresorg(air[,2],air[,1],5)
#
runmeans<-function (x,y,k) 
{
#************************************************************
#Funcion que calcula el suavizador por running means
#Inputs: la variable de respuesta y, la variable predictora x y
#        El numero de vecinos a cada lado de una observacion
#Output: El plot del suavizador superpuesto en el diagrama de dispersion 
#Edgar Acuna, mayo 2003
#*********************************************************
n<-length(x)
rm<-rep(0,n)
for(i in 1:n)
{ind1<-max(i-k,1)
ind2<-min(i+k,n)
tempo<-y[ind1:ind2]
rm[i]<-mean(tempo)
}
plot(x,y)
lines(sort(x),rm,type="l")
title("Running means")
}
#
#Aplicando running means al conjunto de datos air
runmeans(air[,2],air[,1],3)
#
# Regression por kernels
plot(air[,2],air[,1])
lines(ksmooth(air[,2],air[,1],kernel="normal",bandwidth=5))
title("regresion usando kernels")
#
# Lowess regression
plot(air[,2],air[,1])
lines(lowess(air[,2],air[,1]))
title("lowess")