set.seed(127)
err<-NULL
x<-seq(1,200,1)
y<-3.8*x+25+err
for(i in c(1:200)) {
err[i]<-rnorm(1,0,i)
}
y<-3.8*x+25+err;
errno<-NULL
for(i in c(1:200)) {
errno[i]<-rnorm(1,0,15)
};
yn<-3.8*x+25+errno;
naivelm<-lm(y~x)
sqresid<-naivelm$residuals^2
loc<-numeric()
for (i in c(2:200)) {
if (sqresid[i] > max(sqresid[1:i-1]))
loc[i]<-i
}
sqrno<-lm(yn~x)$residuals^2
locx<-numeric()
for (i in c(2:200)) {
if (sqrno[i] > max(sqrno[1:i-1]))
locx[i]<-i
}
#Non para plot
par(mfrow=c(2,1))
plot(x,sqresid/sd(sqresid),type="l",col="blue",ylab="Studentized Residuals",main="Peaks in Goldfeld-Quandt Non-parametric test",sub="Heteroskedastic Errors")
points(loc[!is.na(loc)],sqresid[loc[!is.na(loc)]]/sd(sqresid),pch=20,col="blue")
plot(x,sqrno/sd(sqrno),type="l",col="green",ylab="Studentized Residuals",sub="Homoskedastic Errors")
points(locx[!is.na(locx)],sqrno[locx[!is.na(locx)]]/sd(sqrno),pch=20,col="green")