## A function to graphically examine the residuals

## with a view to looking for and modelling heteroscedasticity.

## This starts off the function. It will be called "heterosc"

heterosc <- function(lmfit)

{

# Check to make sure we are operating on the right sort of object

if(class(lmfit) != "lm")

stop("Data must be from lm")

# Do the studentized residuals calculations.

thii <- lm.influence(lmfit)$hat

ei <- residuals(lmfit)

s <- sqrt(var(residuals(lmfit)))

eip <- ei/(s * sqrt(1 - thii))

studres <- eip * sqrt((lmfit$df - 1)/(lmfit$df - thii^2))

# Now start up the plots.

# This will give a2 by 2 grid of plots

par(mfrow=c(2,2))

plot(predict(lmfit),studres,xlab="Predicted",ylab="Studentized")

# The "lines' command adds lines to an existing plot.

lines(smooth.spline(predict(lmfit),studres),col=3,lwd=5)

plot(predict(lmfit),sqrt(studres^2),

xlab="Predicted",ylab="abs(Studentized)")

lines(smooth.spline(predict(lmfit),sqrt(studres^2)),col=3,lwd=5)

plot(predict(lmfit),studres^2,xlab="Predicted",

ylab="squared(Studentized)")

lines(smooth.spline(predict(lmfit),studres^2),col=3,lwd=5)

qqplot(

qchisq((1:length(predict(lmfit)))/

(length(predict(lmfit))), 1),

studres^2,

ylab="Observed quantile",

xlab = "Expected quantiles from Chisq(1)" )

# This adds thye 45 degree line to the QQ-plot.

abline(c(0,1),col=3,lwd=2)

par(mfrow = c(1,1))

invisible()

}

# Now throw this lot into the gui.

create.menu.hetcheck<- function() {

guiCreate("MenuItem",

Name="SPlusMenuBar$Zinn$hetcheck" ,

Type="MenuItem",

Action="Function",

Command="heterosc",

MenuItemText="&Heteroscedasticity Diagnostics")

invisible()

}

create.menu.hetcheck()

create.menu.hetcheck()

create.toolbar.het <- function(){

guiCreate("Toolbar",Name = "Zinn")

guiCreate("ToolbarButton",

Name = "Zinn$hetcheck",

Type = "Button",

Action = "Function",

Command = "heterosc",

TipText = "Heteroscedasticity check", ShowDialogOnRun=T)

invisible()

}

create.toolbar.het()

#guiRemove("MenuItem",Name="SPlusMenuBar$Zinn$hetcheck")

#guiRemove("MenuItem",Name="SPlusMenuBar$Zinn")

#guiRemove("ToolbarButton",Name="Zinn$hetcheck")

#guiRemove("Toolbar",Name="Zinn")