#
# MATH 3330 Section B   F01

# Week 8

# Chapter 7: regression with 'polytomous' variables
# Uses the 'prestige' data set from Lesson 2

summary(prestige)
xyplot( Income ~ Education | Type, prestige,  panel = panel.xyell)

contrasts( prestige$Type )

fit1 <- lm(Income ~ Education , prestige, na.action = na.omit)
summary(fit1)
anova(fit1)
fit1$contrasts


pred <- expand.grid( Education = 6:16, Type = levels( prestige$Type) )

pred$Inc1 <- predict( fit1, newdata = pred )

xyplot( Inc1 ~ Education , pred , type = 'l')


fit2 <- lm(Income ~ Education + Type, prestige, na.action = na.omit)
summary(fit2)
anova(fit2)

pred$Inc2 <- predict( fit2, newdata = pred )

xyplot( Inc2 ~ Education | Type, pred,  type = 'l')
xyplot( Inc2 ~ Education, pred, groups = Type, panel = panel.superpose, type = 'l')


fit3 <- lm(Income ~ Education * Type, prestige, na.action = na.omit)
summary(fit3)
anova(fit3)

pred$Inc3 <- predict( fit3, newdata = pred )

xyplot( Inc3 ~ Education | Type, pred,  type = 'l')
xyplot( Inc3 ~ Education, pred, groups = Type, panel = panel.superpose, type = 'l')







# example of traditional graphics:

prestige				# data.frame
names(prestige)		# names of variables in data frame
prestige$Income		# a particular variable

# where variables are:

search()				# list of all directories
objects()				# contents of first directory
objects(where=2)		# contents of second directory

# turning a data frame into a directory:

attach(prestige)
search()
objects(where=2)
detach("prestige")		# note that you need quotes here
search()


# example of traditional graphics

plot( prestige$Education, prestige$Income )	# plot( x, y)
text( prestige$Education, prestige$Income, row.names(prestige) )

attach(prestige)
plot(Education, Income)		
identify( Education, Income, row.names(prestige) )

detach("prestige")


# example of Trellis graphics


xyplot( Income ~ Education, prestige)   # plot( y ~ x, data.frame )
identify( xyplot( Income ~ Education , prestige ) )

# a strength of Trellis is easy 'co-plots'

identify( xyplot( Income ~ Education |  Type, prestige ) )

# a weakness is that it's hard to add elements to the plot
# but 'panel' functions give us a way.  It's harder to get
# hang of it at first ...


xyplot(Income ~ Education | Type , prestige, 		
	panel = function(x,y,...) {
		panel.xyplot(x,y, ...)							# plot the data
		panel.lmline(x,y,...)							# plot a least-sq. line
		panel.loess(x, y, span = 2/3,...)			# plot a loess line
	}
)


z <- xyplot(Income ~ Education | Type , prestige, 		
	panel = function(x,y,...) {
		panel.xyplot(x,y, ...)
		panel.lmline(x,y,...)
		panel.loess(x,y, span = 1/3,...)
	}
)

print(z)
identify(z)

z <- xyplot(Income ~ Education | Type , prestige, 		
	panel = function(x,y,...) {
		panel.xyplot(x,y, ...)
		panel.lmline(x,y,...)
		panel.loess(x,y, span = 1/3, family = "symmetric",...)
	}
)

print(z)
identify(z)

z <- xyplot(Income ~ Education | Type , prestige, 		
	panel = function(x,y,...) {
		panel.xyplot(x,y, ...)
		panel.lmline(x,y,...)
		panel.loess(x,y, span = 1/3, family = 'gaussian',...)
		panel.loess(x,y, span = 1/3, family = 'symmetric',lty=2,...)
	}
)
print(z)
identify(z)

z <- xyplot(Income ~ Education | Type , prestige, 		
	panel = function(x,y,...) {
		panel.xyplot(x,y, ...)
		panel.lmline(x,y,...)
		panel.loess(x,y, span = 1/3, family = 'gaussian',...)
		panel.loess(x,y, span = 1/3, family = 'symmetric',lty=2,...) 
		fit <- supsmu(x,y)
		panel.xyplot(fit$x,fit$y, type = 'l',col=3,lty=3)
		
	}
)
print(z)