#
# MATH 3330 B,  October 2001
# Week 4b
# Basic linear regression and connections with
# variance, covariance, correlation
#

## Variance of a matrix

d.mat <- cbind( Education = prestige$Education, Income = prestige$Income )
dim( d.mat )
d.mat[ 1:10, ]

# Variance-Covariance Matrix ( Variance Matrix )

v.mat <- var(d.mat)
v.mat

# Alternate formulas for Regression objects
# See notes

# Visualizing the Variance Matrix:

#
# Circles and ellipses
#

# a circle

rads <- 2* pi * (0:100) / 100			# 100 angles around a circle in radians
circle <- cbind(cos(rads), sin(rads))	# 100 points around a circle as rows of matrix
circle[ 1:10, ]							# starts at (1,0) and goes counterclockwise
plot(circle)

par(pty='s')
plot(circle, type = 'l')

# A linear transformation of a circle

Amat <- cbind ( c(2,1), c(-2, 3))		# a 2 x 2 matrix, i.e. a lin. trans: R2 to R2

plot( circle %*% t(Amat) , type = 'l' )
plot( circle %*% t(Amat) , type = 'l' , xlim=c(-3,3), ylim = c(-3,3))
lines( circle )		# original circle

# Add a unit square to the circle

thing <- rbind( circle, c(NA,NA), c(0,0),c(1,0),c(1,1),c(0,1), c(0,0))
plot( thing , type = 'l' )
plot( thing %*% t(Amat), type = 'l', xlim = c(-4,4), ylim = c(-4,4) )
lines( thing )

#
#
#

data.ell <-  function(x, y, radius = 1,...) {
	#
	# Draws a data ellipse by transforming and translating a circle
	#
	dmat <- cbind( x, y )
	v <- var( dmat )
	m <- apply( dmat, 2, mean)
	rads <- 2* pi * (0:100) / 100
	circle <- cbind(cos(rads), sin(rads))
	lines ( t(m + radius * t( circle %*% chol(v) )),...)
}

#
# Some examples:
#

plot( prestige$Education, prestige$Income )
data.ell ( prestige$Education, prestige$Income )





pairs(prestige)

pairs(prestige, panel = function(x,y,...) {
	panel.xyplot( x, y, ...)
	dd <- na.omit( cbind(x,y))
	data.ell( dd[,1], dd[,2],...)
})


dd <- data.frame( x=rnorm(100))
dd$y <- dd$x + .6* rnorm(100)
plot(dd$x,dd$y)
data.ell( dd$x, dd$y )
abline(coef(lm(y~x,dd)))
data.ell( dd$x, dd$y , rad = 2)