# This file demonstrates linear discriminant function # theory when using a mixture with known probabilities. # # No data is involved, just the densities. # Only two densities: # Linear discriminant function: common covariance matrix # Quadratic discriminant function: different covariance matrices # # Authors: # H.A. Chipman, 2003 # R.W. Oldford, 2004 # # # The following sets up the grid of points on which the densities will # be evaluated. x1.grid <- rep(seq(0,8,l=40),48) x2.grid <- rep(seq(-5,18,l=48),rep(40,48)) x.grid <- cbind(x1.grid,x2.grid) # A handy function to capture the pathname for the files we will load. web441 <- function(x) {paste('http://www.undergrad.math.uwaterloo.ca/~stat441/R-code/', x, sep='') } # get code for creating densities # (in class we just loaded this file since the # machine wasn't connected to the internet). source(web441('dmvnorm.R')) # First specify the multivariate normals parameters for each # density mu0 <- c(3,3) mu1 <- c(6,8) Sigma0 <- cbind(c(1,2),c(2,20)) # To demonstrate effect on linear discriminant we use Sigma1 <- Sigma0 # Or for quadratic we use: (TRY EACH ONE ON ALL THAT FOLLOWS) Sigma1 <- cbind(c(2,3),c(3,16)) # Create the class 0 density den0 <- dmvnorm(x.grid,mu0,Sigma0) # and the class 1 density # den1 <- dmvnorm(x.grid,mu1,Sigma1) # Show the density assuming an equal mixture # First, overlaid the others contour(unique(x1.grid),unique(x2.grid),matrix(den0,40,48),nlevels=6,col='red') # identify the mean as 0 points(mu0[1],mu0[2], pch="0") # Repeat for the class 1 density par(new=T) contour(unique(x1.grid),unique(x2.grid),matrix(den1,40,48),nlevels=6,col='blue') # Now overlay the mixture par(new=T) contour(unique(x1.grid),unique(x2.grid), matrix(0.5 * (den0 + den1),40,48),nlevels=8,lwd=2, col='brown') # Now have the mixture on its own with the discriminant function overlaid. # contour(unique(x1.grid),unique(x2.grid), matrix(0.5 * (den0 + den1),40,48),nlevels=8,lwd=1, col='brown') par(new=T) contour(unique(x1.grid),unique(x2.grid), matrix(log(den0)-log(den1),40,48),nlevels=1,lwd=3,col='black') # What if the data are more likely to come from the class 0 # than from the class 1? # # If they are twice as likely from class 0 than class 1 # then this could be described as a different mixture. # say as contour(unique(x1.grid),unique(x2.grid), matrix(((2/3) * den0 + (1/3)* den1),40,48),nlevels=10,lwd=1, col='brown') # Overlaying the discriminant function as before par(new=T) contour(unique(x1.grid),unique(x2.grid), matrix(log(den0)-log(den1),40,48),nlevels=1,lwd=3,col='black') # # This can be accommodated by making use of these probabilities # as in par(new=T) contour(unique(x1.grid),unique(x2.grid),levels=c(0), matrix(log(den0)-log(den1)+log(2),40,48),nlevels=1,lwd=2,col='green') # # Similarly, had class 0 been half as likely to occur as class 1 # contour(unique(x1.grid),unique(x2.grid), matrix(((1/3) * den0 + (2/3)* den1),40,48),nlevels=10,lwd=1, col='brown') par(new=T) contour(unique(x1.grid),unique(x2.grid), matrix(log(den0)-log(den1),40,48),nlevels=1,lwd=3,col='black') # # This can be accommodated by making use of these probabilities # as in par(new=T) contour(unique(x1.grid),unique(x2.grid),levels=c(0), matrix(log(den0)-log(den1)+log(1/2),40,48),nlevels=1,lwd=2,col='green') # If they are 10 times more likely from class 0 than class 1 # then this could be described as a different mixture. # say overlaid as follows: par(new=T) contour(unique(x1.grid),unique(x2.grid), matrix(((10/11) * den0 + (1/11)* den1),40,48),nlevels=30,lwd=1, col='brown') # note that you can also modify the prior probabilities, which has the # has the effect of shifting the line one way or the other. We add to # the log of the ratio the log of the ratio of prior probs. For # example, if class 0 is twice as likely as class 1, we do: par(new=T) contour(unique(x1.grid),unique(x2.grid),levels=c(0), matrix(log(den0)-log(den1)+log(2),40,48),nlevels=1,lwd=2,col='green') # Or, if class 0 is half as likely as class 1, we do: par(new=T) contour(unique(x1.grid),unique(x2.grid),levels=c(0), matrix(log(den0)-log(den1)+log(1/2),40,48),nlevels=1,lwd=2,col='orange') # Or, if class 0 is 10 times as likely as class 1, we do: par(new=T) contour(unique(x1.grid),unique(x2.grid),levels=c(0), matrix(log(den0)-log(den1)+log(10),40,48),nlevels=1,lwd=2,col='blue')