Difference between revisions of "Applications of Latent Variable Methods"

From Latent Variable Methods in Engineering
Jump to navigation Jump to search
 
(15 intermediate revisions by the same user not shown)
Line 1: Line 1:
== Class 7: soft sensors (21 October)  ==
{{ClassSidebarYouTube
| date = Various dates
| vimeoID1 = PrP_az4viJs
| vimeoID2 = 9Lk3_9A588o
| vimeoID3 =pVlUXzjB1Ec
| vimeoID4 =jzopKpLdEOI
| vimeoID5 =KRJUirDnkdo
| vimeoID6 =yCPHKutbdLs
| vimeoID7 =guk9CYZOHIM
| vimeoID8 =pEyCc1DL4wc
| vimeoID9 =
| course_notes_PDF =
| course_notes_alt = Course notes
| overheads_PDF =
| overheads_PDF_alt = Projector notes
| assignment_instructions =
| assignment_solutions =
| video_download_link_MP4 = http://learnche.mcmaster.ca/media/LVM-2011-Class-07B.mp4
| video_download_link2_MP4 = http://learnche.mcmaster.ca/media/LVM-2011-Class-07C.mp4
| video_download_link3_MP4 = http://learnche.mcmaster.ca/media/LVM-2011-Class-08A.mp4
| video_download_link4_MP4 = http://learnche.mcmaster.ca/media/LVM-2011-Class-08B.mp4
| video_download_link5_MP4 = http://learnche.mcmaster.ca/media/LVM-2011-Class-08C.mp4
| video_download_link6_MP4 = http://learnche.mcmaster.ca/media/LVM-2011-Class-10A.mp4
| video_download_link7_MP4 = http://learnche.mcmaster.ca/media/LVM-2011-Class-10B.mp4
| video_download_link8_MP4 = http://learnche.mcmaster.ca/media/LVM-2011-Class-10C.mp4
| video_download_link9_MP4 =
| video_download_link_MP4_size = 143 Mb
| video_download_link2_MP4_size = 178 Mb
| video_download_link3_MP4_size = 251 Mb
| video_download_link4_MP4_size = 318 Mb
| video_download_link5_MP4_size = 187 Mb
| video_download_link6_MP4_size = 205 Mb
| video_download_link7_MP4_size = 293 Mb
| video_download_link8_MP4_size = 299 Mb
| video_download_link9_MP4_size =  Mb
| video_notes1 =
| video_notes2 =
| video_notes3 =
| video_notes4 =
| video_notes5 =
| video_notes6 =
}}
== Class 7: soft sensors (21 October 2011)  ==
 
[[Image:Nuvola_mimetypes_pdf.png|20px|link=Media:Lvm-class-7.pdf]] [[Media:Lvm-class-7.pdf|Download the class slides]] (PDF)


<pdfreflow>   
class_date        = 21 October 2011  [3.6 Mb]
button_label      = Create my projector slides!
show_page_layout  = 1
show_frame_option = 1
pdf_file          = lvm-class-7.pdf
</pdfreflow>


Please download and bring the following data sets to class:
Please download and bring the following data sets to class:
* Sawdust: http://datasets.connectmv.com/info/sawdust
* Sawdust: http://openmv.net/info/sawdust
* Distillation tower: http://datasets.connectmv.com/info/distillation-tower
* Distillation tower: http://openmv.net/info/distillation-tower
 


== Class 8: advanced preprocessing, empirical models, adaptive models (21 October)  ==
== Class 8: advanced preprocessing, empirical models, adaptive models (21 October 2011)  ==


<pdfreflow>   
[[Image:Nuvola_mimetypes_pdf.png|20px|link=Media:Lvm-class-8.pdf]] [[Media:Lvm-class-8.pdf|Download the class slides]] (PDF)
class_date        = 28 October 2011  [1.3 Mb]
button_label      = Create my projector slides!
show_page_layout  = 1
show_frame_option = 1
pdf_file          = lvm-class-8.pdf
</pdfreflow>


Please download and bring the following data sets to class:
Please download and bring the following data sets to class:
* Kamyr digester: http://datasets.connectmv.com/info/kamyr-digester
* Kamyr digester: http://openmv.net/info/kamyr-digester


'''Note''': the first few slides are similar to last week's slides that we didn't cover in class.
'''Note''': the first few slides are similar to last week's slides that we didn't cover in class.
== Class 9: Dealing with image data (04 November 2011)  ==
See the [[Dealing_with_image_data|Dealing with image data]] section.
== Class 10: Classification (11 November 2011)  ==
[[Image:Nuvola_mimetypes_pdf.png|20px|link=Media:Lvm-class-10.pdf]] [[Media:Lvm-class-10.pdf|Download the class slides]] (PDF)
Also download and bring the following data set to class:
* [http://learnche.org/images/a/ae/Olive-oil.csv Olive oil]


=== R script used to compare models ===
=== R script used to compare models ===
<syntaxhighlight lang="R">
<syntaxhighlight lang="R">
# Load the data and create training/testing splits
# Load the data and create training/testing splits
cheese <- read.csv('http://datasets.connectmv.com/file/cheddar-cheese.csv')
cheese <- read.csv('http://openmv.net/file/cheddar-cheese.csv')
N = dim(cheese)[1]
N = dim(cheese)[1]
set.seed(2)
set.seed(2)
Line 102: Line 143:
# Predictions, using test set swap
# Predictions, using test set swap
model.nn1 <- neuralnet(Taste ~ Acetic + H2S + Lactic + Random, cheese[part1,], hidden=1, linear.output = TRUE)  
model.nn1 <- neuralnet(Taste ~ Acetic + H2S + Lactic + Random, cheese[part1,], hidden=1, linear.output = TRUE)  
# model.nn2 <- neuralnet(Taste ~ Acetic + H2S + Lactic + Random, cheese[part2,], hidden=1, linear.output = TRUE, algorithm="rprop+", rep=3)   
# model.nn2 <- neuralnet(Taste ~ Acetic + H2S + Lactic + Random, cheese[part2,], hidden=1, linear.output = TRUE)   
# Previous model doesn't converge with that particular training set. Add one extra observation to the part 2 data: converges now
# Previous model doesn't converge with that particular training set. Add one extra observation to the part 2 data: converges now
model.nn2 <- neuralnet(Taste ~ Acetic + H2S + Lactic + Random, cheese[c(1, part2),], hidden=1, linear.output = TRUE)   
model.nn2 <- neuralnet(Taste ~ Acetic + H2S + Lactic + Random, cheese[c(1, part2),], hidden=1, linear.output = TRUE)   

Latest revision as of 14:06, 17 September 2018

Class date(s): Various dates
Video material (part 1)
Download video: Link (plays in Google Chrome) [143 Mb]


Video material(part 2)
Download video: Link (plays in Google Chrome) [178 Mb]


Video material (part 3)
Download video: Link (plays in Google Chrome) [251 Mb]


Video material (part 4)
Download video: Link (plays in Google Chrome) [318 Mb]


Video material (part 5)
Download video: Link (plays in Google Chrome) [187 Mb]


Video material (part 6)
Download video: Link (plays in Google Chrome) [205 Mb]


Video material (part 7)
Download video: Link (plays in Google Chrome) [293 Mb]


Video material (part 8)
Download video: Link (plays in Google Chrome) [299 Mb]

Class 7: soft sensors (21 October 2011)

Nuvola mimetypes pdf.png Download the class slides (PDF)


Please download and bring the following data sets to class:

Class 8: advanced preprocessing, empirical models, adaptive models (21 October 2011)

Nuvola mimetypes pdf.png Download the class slides (PDF)

Please download and bring the following data sets to class:

Note: the first few slides are similar to last week's slides that we didn't cover in class.

Class 9: Dealing with image data (04 November 2011)

See the Dealing with image data section.

Class 10: Classification (11 November 2011)

Nuvola mimetypes pdf.png Download the class slides (PDF)

Also download and bring the following data set to class:

R script used to compare models

# Load the data and create training/testing splits
cheese <- read.csv('http://openmv.net/file/cheddar-cheese.csv')
N = dim(cheese)[1]
set.seed(2)
cheese$Random <- rnorm(N, 1)
part1 <- seq(1, dim(cheese)[1], 2)
part2 <- seq(2, dim(cheese)[1], 2)
cols <- c("Acetic", "H2S", "Lactic", "Random")

# Correlation between variables
cor(cheese)

# Create a new cheese with these properties
new.cheese <- data.frame(Acetic=9.0, H2S=4.0, Lactic=1.8, Random=0)
write.csv(cheese, 'cheese-with-random-data.csv')  # to double check calculations in other software packages

# Scatter plot matrix
# -----------------------
library(car)
bitmap('cheese-plots.png', type="png256", width=6, height=6, res=300, pointsize=14)
par(mar=c(1.5, 1.5, 1.5, 0.5))  # (bottom, left, top, right); defaults are par(mar=c(5, 4, 4, 2) + 0.1)
par(cex.lab=1.2, cex.main=1.2, cex.sub=1.2, cex.axis=1.2)
scatterplotMatrix(cheese[,2:6], col=c(1,1,1), smooth=FALSE)
dev.off()

# Multiple linear regression model
# -----------------------
model.lm <- lm(Taste ~ Acetic + H2S + Lactic + Random, data=cheese)
summary(model.lm)
confint(model.lm)
predict(model.lm, new.cheese)  # 29.9
RMSEE.lm <- sqrt(mean((cheese$Taste - predict(model.lm))**2))
# Predictions using the "test set swap" procedure
model.lm1 <- lm(Taste ~ Acetic + H2S + Lactic, data=cheese[part1,])
model.lm2 <- lm(Taste ~ Acetic + H2S + Lactic, data=cheese[part2,])
RMSEP.lm1 <- sqrt(mean((cheese$Taste[part1] - predict(model.lm2, cheese[part1,]))**2))
RMSEP.lm2 <- sqrt(mean((cheese$Taste[part2] - predict(model.lm1, cheese[part2,]))**2))
RMSEP.lm <- mean(c(RMSEP.lm1, RMSEP.lm2))

# Ordinary Least Squares model: H2S selected since it has highest correlation with Taste
# ----------------------------
model.lm.H2S <- lm(Taste ~ H2S, data=cheese)
summary(model.lm.H2S)
confint(model.lm.H2S)
predict(model.lm.H2S, new.cheese) 
RMSEE.lm.H2S <- sqrt(mean((cheese$Taste - predict(model.lm.H2S))**2))
# Predictions
model.H2S.lm1 <- lm(Taste ~ H2S, data=cheese[part1,])
model.H2S.lm2 <- lm(Taste ~ H2S, data=cheese[part2,])
RMSEP.H2S.lm1 <- sqrt(mean((cheese$Taste[part1] - predict(model.H2S.lm2, cheese[part1,]))**2))
RMSEP.H2S.lm2 <- sqrt(mean((cheese$Taste[part2] - predict(model.H2S.lm1, cheese[part2,]))**2))
RMSEP.H2S.lm <- mean(c(RMSEP.H2S.lm1, RMSEP.H2S.lm2))


#Neural network model
#--------------------
library(neuralnet)
model.nn <- neuralnet(Taste ~ Acetic + H2S + Lactic + Random, cheese, 
                      hidden=1, 
                      likelihood=TRUE,
                      algorithm="rprop+",
                      linear.output = TRUE)  # because we want predictions of y, not a 0/1 classifier
yhat <- compute(model.nn, cheese[,c(2,3,4,6)])$net.result 
yhat <- model.nn$net.result[[1]][,1]  # two ways of getting the neural network predictions
resid <- cheese$Taste - yhat
RMSEE.nn <- sqrt(mean(resid**2))
R2.nn <- 1 - var(resid) / var(cheese$Taste )
compute(model.nn, new.cheese)   # 38.8

# Predictions, using test set swap
model.nn1 <- neuralnet(Taste ~ Acetic + H2S + Lactic + Random, cheese[part1,], hidden=1, linear.output = TRUE) 
# model.nn2 <- neuralnet(Taste ~ Acetic + H2S + Lactic + Random, cheese[part2,], hidden=1, linear.output = TRUE)   
# Previous model doesn't converge with that particular training set. Add one extra observation to the part 2 data: converges now
model.nn2 <- neuralnet(Taste ~ Acetic + H2S + Lactic + Random, cheese[c(1, part2),], hidden=1, linear.output = TRUE)   
yhat.nn.2 <- compute(model.nn1, cheese[part2, c(2,3,4,6)])$net.result   # predict part2's data, given model built from part 1
yhat.nn.1 <- compute(model.nn2, cheese[part1, c(2,3,4,6)])$net.result   # predict part1's data, given model built from part 2
RMSEP.nn.2 <- sqrt(mean((cheese$Taste[part2] - yhat.nn.2)**2))  # 14.9
RMSEP.nn.1 <- sqrt(mean((cheese$Taste[part1] - yhat.nn.1)**2))
RMSEP.nn <- mean(c(RMSEP.nn.1, RMSEP.nn.2))

# PCR model (one component): long way (manually create PCR model, then linear regression model)
# -----------------------
model.pca <- prcomp(cheese[,cols], scale=TRUE)
summary(model.pca)
P <- model.pca$rotation
T <- model.pca$x
model.pcr <- lm(cheese$Taste ~ T[,1])
summary(model.pcr)
RMSEE.pcr <- sqrt(mean((cheese$Taste - predict(model.pcr))**2))

# regression coefficients
P[,1] * model.pcr$coef[2] / model.pca$scale

# Model for part 1
model.pca.1 <- prcomp(cheese[part1,cols], scale=TRUE)
P.1 <- model.pca.1$rotation
T.1 <- model.pca.1$x
model.pcr.1 <- lm(cheese$Taste[part1] ~ T.1[,1])
# Model for part 2
model.pca.2 <- prcomp(cheese[part2,cols], scale=TRUE)
P.2 <- model.pca.2$rotation
T.2 <- model.pca.2$x
model.pcr.2 <- lm(cheese$Taste[part2] ~ T.2[,1])

# Center and scale the new data, using the training set's centering and scaling vectors
x.2.new <- scale(cheese[part2,cols], center= model.pca.1$center, scale=model.pca.1$scale)
x.1.new <- scale(cheese[part1,cols], center= model.pca.2$center, scale=model.pca.2$scale)
# Project the testing data (.new matrices) onto the opposite model
T.2.hat <- x.2.new %*% P.1
T.1.hat <- x.1.new %*% P.2

predict(model.pcr.1)
# Predictions
pred.2 <- model.pcr.1$coef[1] + model.pcr.1$coef[2] * T.2.hat[,1]
pred.1 <- model.pcr.2$coef[1] + model.pcr.2$coef[2] * T.1.hat[,1]
RMSEP.pcr.lm1 <- sqrt(mean((cheese$Taste[part1] - pred.1)**2))
RMSEP.pcr.lm2 <- sqrt(mean((cheese$Taste[part2] - pred.2)**2))
RMSEP.pcr <- mean(c(RMSEP.pcr.lm1, RMSEP.pcr.lm2))

# PCR using the "pls" package (gets the same results as the longer code above)
# ----------------------------
library(pls)
model.pcr.quick <- mvr(Taste ~ Acetic + H2S + Lactic + Random, data=cheese, method="svdpc", scale=TRUE)
summary(model.pcr.quick)
# R2_taste, using 1PC = 62.04%
RMSEE.pcr <- sqrt(mean((cheese$Taste - predict(model.pcr.quick)[,,1])**2))
predict.mvr(model.pcr.quick, new.cheese)  # 54.3


model.pcr.quick.1 <- mvr(Taste ~ Acetic + H2S + Lactic + Random, data=cheese, method="svdpc", scale=TRUE, subset=part1)
model.pcr.quick.2 <- mvr(Taste ~ Acetic + H2S + Lactic + Random, data=cheese, method="svdpc", scale=TRUE, subset=part2)

RMSEP.pcr.quick1 <- sqrt(mean((cheese$Taste[part1] - predict(model.pcr.quick.2, cheese[part1,])[,,1])**2))
RMSEP.pcr.quick2 <- sqrt(mean((cheese$Taste[part2] - predict(model.pcr.quick.1, cheese[part2,])[,,1])**2))
RMSEP.pcr <- mean(c(RMSEP.pcr.quick1, RMSEP.pcr.quick2))

# Regression coefficients with 1 PC
coef(model.pcr.quick, ncomp=1)

# PLS using the "pls" package
# ----------------------------
library(pls)
model.pls <- mvr(Taste ~ Acetic + H2S + Lactic + Random, data=cheese, method="simpls", scale=TRUE)
summary(model.pls)
# R2_taste, using 1PC = 62.7%
RMSEE.pls <- sqrt(mean((cheese$Taste - predict(model.pls)[,,1])**2))
predict.mvr(model.pls, new.cheese) 

model.pls.1 <- mvr(Taste ~ Acetic + H2S + Lactic + Random, data=cheese, method="simpls", scale=TRUE, subset=part1)
model.pls.2 <- mvr(Taste ~ Acetic + H2S + Lactic + Random, data=cheese, method="simpls", scale=TRUE, subset=part2)

RMSEP.pls.1 <- sqrt(mean((cheese$Taste[part1] - predict(model.pls.2, cheese[part1,])[,,1])**2))
RMSEP.pls.2 <- sqrt(mean((cheese$Taste[part2] - predict(model.pls.1, cheese[part2,])[,,1])**2))
RMSEP.pls <- mean(c(RMSEP.pls.1, RMSEP.pls.2))

# Coefficients and reliability intervals from ProMV 11.08,  but divide them by "model.pls$scale" and times by sd(y)
#               Coeff               Lower bound         Upper bound
# H2S       1   0.341649600476995   0.214433534821241   0.46886566613275
# Lactic    2   0.318015562701819   0.237906914060081   0.398124211343557
# Acetic    3   0.248419848376738   0.176929971367371   0.319909725386105
# Random	4	0.126124729088615	-0.034083388621767	0.286332846798997

c(0.3416, 0.3180, 0.2484, 0.1261) / model.pls$scale * sd(cheese$Taste)