Applied Exercises of Chapter 4
Question 10
10a
library(ISLR)
str(Weekly)'data.frame': 1089 obs. of 9 variables:
$ Year : num 1990 1990 1990 1990 1990 1990 1990 1990 1990 1990 ...
$ Lag1 : num 0.816 -0.27 -2.576 3.514 0.712 ...
$ Lag2 : num 1.572 0.816 -0.27 -2.576 3.514 ...
$ Lag3 : num -3.936 1.572 0.816 -0.27 -2.576 ...
$ Lag4 : num -0.229 -3.936 1.572 0.816 -0.27 ...
$ Lag5 : num -3.484 -0.229 -3.936 1.572 0.816 ...
$ Volume : num 0.155 0.149 0.16 0.162 0.154 ...
$ Today : num -0.27 -2.576 3.514 0.712 1.178 ...
$ Direction: Factor w/ 2 levels "Down","Up": 1 1 2 2 2 1 2 2 2 1 ...summary(Weekly) Year Lag1 Lag2 Lag3 Lag4
Min. :1990 Min. :-18.1950 Min. :-18.1950 Min. :-18.1950 Min. :-18.1950
1st Qu.:1995 1st Qu.: -1.1540 1st Qu.: -1.1540 1st Qu.: -1.1580 1st Qu.: -1.1580
Median :2000 Median : 0.2410 Median : 0.2410 Median : 0.2410 Median : 0.2380
Mean :2000 Mean : 0.1506 Mean : 0.1511 Mean : 0.1472 Mean : 0.1458
3rd Qu.:2005 3rd Qu.: 1.4050 3rd Qu.: 1.4090 3rd Qu.: 1.4090 3rd Qu.: 1.4090
Max. :2010 Max. : 12.0260 Max. : 12.0260 Max. : 12.0260 Max. : 12.0260
Lag5 Volume Today Direction
Min. :-18.1950 Min. :0.08747 Min. :-18.1950 Down:484
1st Qu.: -1.1660 1st Qu.:0.33202 1st Qu.: -1.1540 Up :605
Median : 0.2340 Median :1.00268 Median : 0.2410
Mean : 0.1399 Mean :1.57462 Mean : 0.1499
3rd Qu.: 1.4050 3rd Qu.:2.05373 3rd Qu.: 1.4050
Max. : 12.0260 Max. :9.32821 Max. : 12.0260 plot(Weekly)
The Volume increased with the Year. No other patterns are clear enough to be observed.
10b
lrf <- glm(Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5 + Volume, data = Weekly, family = 'binomial')
summary(lrf)
Call:
glm(formula = Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5 +
Volume, family = "binomial", data = Weekly)
Deviance Residuals:
Min 1Q Median 3Q Max
-1.6949 -1.2565 0.9913 1.0849 1.4579
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 0.26686 0.08593 3.106 0.0019 **
Lag1 -0.04127 0.02641 -1.563 0.1181
Lag2 0.05844 0.02686 2.175 0.0296 *
Lag3 -0.01606 0.02666 -0.602 0.5469
Lag4 -0.02779 0.02646 -1.050 0.2937
Lag5 -0.01447 0.02638 -0.549 0.5833
Volume -0.02274 0.03690 -0.616 0.5377
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 1496.2 on 1088 degrees of freedom
Residual deviance: 1486.4 on 1082 degrees of freedom
AIC: 1500.4
Number of Fisher Scoring iterations: 4Lag2 is statistically significant, with p-value 3%.
10c
lrp <- predict(lrf, type = 'response')
lrpred <- rep('Down', length(lrp))
lrpred[lrp > 0.5] <- 'Up'
table(lrpred, Weekly$Direction)
lrpred Down Up
Down 54 48
Up 430 557The overall error rate: \((48 + 430) \div (54 + 48 + 430 + 557) = 0.439\)
Logistic regression predict correctly most of time when direction is “Up”: \(557 \div (557 + 48)\), wrongly most of time when it’s “Down”: \(54 \div (54 + 430)\).
10d
isTest <- Weekly$Year > 2008
training <- Weekly[!isTest, ]
testing <- Weekly[isTest, ]
f10d <- glm(Direction ~ Lag2, data = training, family = 'binomial')
pred10d <- predict(f10d, testing, type = 'response')
res10d <- rep('Down', length(pred10d))
res10d[pred10d > 0.5] = 'Up'
table(res10d, testing$Direction)
res10d Down Up
Down 9 5
Up 34 56mean(res10d == testing$Direction)[1] 0.62510e
library(MASS)
f10e <- lda(Direction ~ Lag2, data = training)
pred10e <- predict(f10e, testing)
table(pred10e$class, testing$Direction)
Down Up
Down 9 5
Up 34 56mean(pred10e$class == testing$Direction)[1] 0.62510f
f10f <- qda(Direction ~ Lag2, data = training)
pred10f <- predict(f10f, testing)
table(pred10f$class, testing$Direction)
Down Up
Down 0 0
Up 43 61mean(pred10f$class == testing$Direction)[1] 0.586538510g
library(class)
set.seed(1)
f10g <- knn(as.matrix(training$Lag2), as.matrix(testing$Lag2), training$Direction, k = 1)
table(f10g, testing$Direction)
f10g Down Up
Down 21 30
Up 22 31mean(f10g == testing$Direction)[1] 0.5Notice the first 2 parameters of knn() function must be data frame or matrix.
对本章第8题的验证： KNN 算法在 \(K=1\) 时的训练正确率是否为 100%？
library(tidyverse)
set.seed(1)
training <- Weekly[!isTest, ]
mk1 <- knn(as.matrix(training$Lag2), as.matrix(training$Lag2), training$Direction, k = 1)
table(mk1, training$Direction)
mk1 Down Up
Down 422 22
Up 19 522mean(mk1 == training$Direction)[1] 0.9583756正确率居然不是 100%？
为什么会有41个点的 Direction 值与自己不同？ 让我们先把它们找出来：
outliers <- which(mk1 != training$Direction)
length(outliers)[1] 41outliers[1:2][1] 12 108异常点共有41个，与上面 table() 的输出 22 + 19 一致。 下面分析一下这些异常点和正常点有什么不一样。
首先取前两个异常点（下标分别为12和108）作为分析对象， 取它们后面的值作为对比，观察它们各自在半径0.01范围内的邻居：
t1 <- 12
t2 <- 108
eps <- 0.01
filter(training, Lag2 > training$Lag2[t1] - eps & Lag2 < training$Lag2[t1] + eps)filter(training, Lag2 > training$Lag2[t1 + 1] - eps & Lag2 < training$Lag2[t1 + 1] + eps)filter(training, Lag2 > training$Lag2[t2] - eps & Lag2 < training$Lag2[t2] + eps)filter(training, Lag2 > training$Lag2[t2 + 1] - eps & Lag2 < training$Lag2[t2 + 1] + eps)不难看出，这些点之所以与“自己”的 Direction 取值不同，是因为有一个或多个相同值的点存在导致的。
最后回答第8题引出的问题：
当数据集里没有重叠点时，\(K=1\)的
knn在自身训练集上的正确率是 100%，如果数据集里有重叠点，是否出出现异常点取决于重叠点的标记（这里是 Direction）取值是否一致， 异常点越多，正确率越低。
10h
Logistic Regression in (d) and LDA method in (e) provide the best result (correct rate: 62.5%) on this data.
Question 11
11a
auto11 <- Auto
auto11$mpg01 <- 0
auto11$mpg01[Auto$mpg > median(Auto$mpg)] <- 111b
plot(auto11)
par(mfrow=c(3,3))
boxplot(cylinders ~ mpg01, data = auto11)
boxplot(displacement ~ mpg01, data = auto11)
boxplot(horsepower ~ mpg01, data = auto11)
boxplot(weight ~ mpg01, data = auto11)
boxplot(acceleration ~ mpg01, data = auto11)
boxplot(year ~ mpg01, data = auto11)
boxplot(origin ~ mpg01, data = auto11)
par(mfrow = c(1,1))
cor(auto11[, -9]) # calculate cor matrix except *name* columncylinder, displacement, horsepower and weight have relative high negative corelations with mpg01 (absolute value are all above 0.65).
11c
isTraining <- sample(nrow(auto11), nrow(auto11) * 0.8)
training <- auto11[isTraining, ]
testing <- auto11[-isTraining, ]11d
f11d <- lda(mpg01 ~ cylinders + displacement + horsepower + weight, data = auto11)
pred11d <- predict(f11d, testing)
mean(pred11d$class == testing$mpg01)The test error is \(1 - 0.873 = 0.127%\).
11e
f11e <- qda(mpg01 ~ cylinders + displacement + horsepower + weight, data = auto11)
pred11e <- predict(f11e, testing)
1 - mean(pred11e$class == testing$mpg01)11f
f11f <- glm(mpg01 ~ cylinders + displacement + horsepower + weight, data = auto11)
pred11f <- predict(f11f, testing)
res11f <- rep(0, length(pred11f))
res11f[pred11f > 0.5] = 1
1 - mean(res11f == testing$mpg01)11g
set.seed(1)
train11g <- auto11[isTraining, c('cylinders', 'displacement', 'horsepower', 'weight')]
test11g <- auto11[-isTraining, c('cylinders', 'displacement', 'horsepower', 'weight')]
trainRes <- auto11$mpg01[isTraining]
f11g1 <- knn(train11g, test11g, trainRes, k = 1)
mean(f11g1 == auto11$mpg01[-isTraining])
f11g3 <- knn(train11g, test11g, trainRes, k = 3)
mean(f11g3 == auto11$mpg01[-isTraining])
f11g5 <- knn(train11g, test11g, trainRes, k = 5)
mean(f11g5 == auto11$mpg01[-isTraining])
f11g100 <- knn(train11g, test11g, trainRes, k = 100)
mean(f11g100 == auto11$mpg01[-isTraining])The minimum error rate is \(1 - 0.937 = 0.063\). \(k = 5\) performs best on this data set.
Question 12
12a
Power <- function() { 2 ^ 3 }
Power()12b & 12c
Power2 <- function(x, a) { x ^ a }
Power2(3, 8)
Power2(10, 3)
Power2(8, 17)
Power2(131, 3)12d
As the current version of R, the last expression is the return value of the function. So it’s unnecessary to use return() function explicitly:
version12e
par(mfrow=c(1,2))
x <- 1:10
plot(x, Power2(x, 2), main = 'f(x) = x^2', xlab = 'x = 1:10', ylab = 'x^2')
plot(log(x), log(Power2(x, 2)), main = 'log(f(x))', xlab = 'x = 1:10', ylab = '2 * log(x)')12f
PlotPower <- function(x, p) {
plot(x, x ^ p, main = 'power p of x')
}
PlotPower(1:10, 3)Question 13
Extract training and testing data, and predict based on the first group of predictors with logistic regression method:
boston13 <- Boston
boston13$crime.rate <- 0
boston13$crime.rate[Boston$crim > median(Boston$crim)] <- 1
isTraining <- sample(nrow(boston13), 0.8 * nrow(boston13))
training <- boston13[isTraining, ]
testing <- boston13[-isTraining, ]
lr.fit1 <- glm(crime.rate ~ zn + indus + chas + nox, data = boston13, family = 'binomial')
lr.pred1 <- predict(lr.fit1, testing)
lr.res1 <- rep(0, length(lr.pred1))
lr.res1[lr.pred1 > median(Boston$crim)] <- 1
mean(lr.res1 == testing$crime.rate)Using another group of predictors:
lr.fit2 <- glm(crime.rate ~ rm + age + dis + rad, data = boston13, family = 'binomial')
lr.pred2 <- predict(lr.fit2, testing)
lr.res2 <- rep(0, length(lr.pred2))
lr.res2[lr.pred2 > median(Boston$crim)] <- 1
mean(lr.res2 == testing$crime.rate)Predict on the 1st group predictors with LDA:
lda.fit1 <- lda(crime.rate ~ zn + indus + chas + nox, data = boston13)
lda.pred1 <- predict(lda.fit1, testing)
mean(lda.pred1$class == testing$crime.rate)Predict on the 2nd group predictors with LDA:
lda.fit2 <- lda(crime.rate ~ rm + age + dis + rad, data = boston13)
lda.pred2 <- predict(lda.fit2, testing)
mean(lda.pred2$class == testing$crime.rate)Predict on the 1st group predictors with KNN (k = 5):
ktrain1 <- boston13[isTraining, c('zn', 'indus', 'chas', 'nox')]
ktest1 <- boston13[-isTraining, c('zn', 'indus', 'chas', 'nox')]
ktrain.result <- boston13$crime.rate[isTraining]
kmod1 <- knn(ktrain1, ktest1, ktrain.result, k = 5)
mean(kmod1 == boston13$crime.rate[-isTraining])Predict on the 2nd group predictors with KNN (k = 5):
ktrain2 <- boston13[isTraining, c('rm', 'age', 'dis', 'rad')]
ktest2 <- boston13[-isTraining, c('rm', 'age', 'dis', 'rad')]
kmod2 <- knn(ktrain2, ktest2, ktrain.result, k = 5)
mean(kmod2 == boston13$crime.rate[-isTraining])KNN (k=5) based on the first group of predictors has the highest accuracy.
LS0tCnRpdGxlOiAiQXBwbGllZCBFeGVyY2lzZXMgb2YgQ2hhcHRlciA0IgpvdXRwdXQ6IGh0bWxfbm90ZWJvb2sKLS0tCgojIFF1ZXN0aW9uIDEwCgojIyAxMGEKCmBgYHtyfQpsaWJyYXJ5KElTTFIpCnN0cihXZWVrbHkpCnN1bW1hcnkoV2Vla2x5KQpwbG90KFdlZWtseSkKYGBgClRoZSAqVm9sdW1lKiBpbmNyZWFzZWQgd2l0aCB0aGUgKlllYXIqLgpObyBvdGhlciBwYXR0ZXJucyBhcmUgY2xlYXIgZW5vdWdoIHRvIGJlIG9ic2VydmVkLgoKIyMgMTBiCgpgYGB7cn0KbHJmIDwtIGdsbShEaXJlY3Rpb24gfiBMYWcxICsgTGFnMiArIExhZzMgKyBMYWc0ICsgTGFnNSArIFZvbHVtZSwgZGF0YSA9IFdlZWtseSwgZmFtaWx5ID0gJ2Jpbm9taWFsJykKc3VtbWFyeShscmYpCmBgYAoqTGFnMiogaXMgc3RhdGlzdGljYWxseSBzaWduaWZpY2FudCwgd2l0aCBwLXZhbHVlIDMlLgoKIyMgMTBjCgpgYGB7cn0KbHJwIDwtIHByZWRpY3QobHJmLCB0eXBlID0gJ3Jlc3BvbnNlJykKbHJwcmVkIDwtIHJlcCgnRG93bicsIGxlbmd0aChscnApKQpscnByZWRbbHJwID4gMC41XSA8LSAnVXAnCnRhYmxlKGxycHJlZCwgV2Vla2x5JERpcmVjdGlvbikKYGBgCgpUaGUgb3ZlcmFsbCBlcnJvciByYXRlOiAkKDQ4ICsgNDMwKSBcZGl2ICg1NCArIDQ4ICsgNDMwICsgNTU3KSA9IDAuNDM5JAoKTG9naXN0aWMgcmVncmVzc2lvbiBwcmVkaWN0IGNvcnJlY3RseSBtb3N0IG9mIHRpbWUgd2hlbiAqZGlyZWN0aW9uKiBpcyAiVXAiOiAkNTU3IFxkaXYgKDU1NyArIDQ4KSQsIHdyb25nbHkgbW9zdCBvZiB0aW1lIHdoZW4gaXQncyAiRG93biI6ICQ1NCBcZGl2ICg1NCArIDQzMCkkLgoKIyMgMTBkCgpgYGB7cn0KaXNUZXN0IDwtIFdlZWtseSRZZWFyID4gMjAwOAp0cmFpbmluZyA8LSBXZWVrbHlbIWlzVGVzdCwgXQp0ZXN0aW5nIDwtIFdlZWtseVtpc1Rlc3QsIF0KZjEwZCA8LSBnbG0oRGlyZWN0aW9uIH4gTGFnMiwgZGF0YSA9IHRyYWluaW5nLCBmYW1pbHkgPSAnYmlub21pYWwnKQpwcmVkMTBkIDwtIHByZWRpY3QoZjEwZCwgdGVzdGluZywgdHlwZSA9ICdyZXNwb25zZScpCnJlczEwZCA8LSByZXAoJ0Rvd24nLCBsZW5ndGgocHJlZDEwZCkpCnJlczEwZFtwcmVkMTBkID4gMC41XSA9ICdVcCcKdGFibGUocmVzMTBkLCB0ZXN0aW5nJERpcmVjdGlvbikKbWVhbihyZXMxMGQgPT0gdGVzdGluZyREaXJlY3Rpb24pCmBgYAoKIyMgMTBlCgpgYGB7cn0KbGlicmFyeShNQVNTKQpmMTBlIDwtIGxkYShEaXJlY3Rpb24gfiBMYWcyLCBkYXRhID0gdHJhaW5pbmcpCnByZWQxMGUgPC0gcHJlZGljdChmMTBlLCB0ZXN0aW5nKQp0YWJsZShwcmVkMTBlJGNsYXNzLCB0ZXN0aW5nJERpcmVjdGlvbikKbWVhbihwcmVkMTBlJGNsYXNzID09IHRlc3RpbmckRGlyZWN0aW9uKQpgYGAKCiMjIDEwZgoKYGBge3J9CmYxMGYgPC0gcWRhKERpcmVjdGlvbiB+IExhZzIsIGRhdGEgPSB0cmFpbmluZykKcHJlZDEwZiA8LSBwcmVkaWN0KGYxMGYsIHRlc3RpbmcpCnRhYmxlKHByZWQxMGYkY2xhc3MsIHRlc3RpbmckRGlyZWN0aW9uKQptZWFuKHByZWQxMGYkY2xhc3MgPT0gdGVzdGluZyREaXJlY3Rpb24pCmBgYAoKIyMgMTBnCgpgYGB7cn0KbGlicmFyeShjbGFzcykKc2V0LnNlZWQoMSkKZjEwZyA8LSBrbm4oYXMubWF0cml4KHRyYWluaW5nJExhZzIpLCBhcy5tYXRyaXgodGVzdGluZyRMYWcyKSwgdHJhaW5pbmckRGlyZWN0aW9uLCBrID0gMSkKdGFibGUoZjEwZywgdGVzdGluZyREaXJlY3Rpb24pCm1lYW4oZjEwZyA9PSB0ZXN0aW5nJERpcmVjdGlvbikKYGBgCgpOb3RpY2UgdGhlIGZpcnN0IDIgcGFyYW1ldGVycyBvZiBga25uKClgIGZ1bmN0aW9uIG11c3QgYmUgZGF0YSBmcmFtZSBvciBtYXRyaXguCgrlr7nmnKznq6DnrKw46aKY55qE6aqM6K+B77yaCktOTiDnrpfms5XlnKggJEs9MSQg5pe255qE6K6t57uD5q2j56Gu546H5piv5ZCm5Li6IDEwMCXvvJ8KYGBge3J9CmxpYnJhcnkodGlkeXZlcnNlKQpzZXQuc2VlZCgxKQp0cmFpbmluZyA8LSBXZWVrbHlbIWlzVGVzdCwgXQptazEgPC0ga25uKGFzLm1hdHJpeCh0cmFpbmluZyRMYWcyKSwgYXMubWF0cml4KHRyYWluaW5nJExhZzIpLCB0cmFpbmluZyREaXJlY3Rpb24sIGsgPSAxKQp0YWJsZShtazEsIHRyYWluaW5nJERpcmVjdGlvbikKbWVhbihtazEgPT0gdHJhaW5pbmckRGlyZWN0aW9uKQpgYGAKCuato+ehrueOh+WxheeEtuS4jeaYryAxMDAl77yfCgrkuLrku4DkuYjkvJrmnIk0MeS4queCueeahCAqRGlyZWN0aW9uKiDlgLzkuI7oh6rlt7HkuI3lkIzvvJ8K6K6p5oiR5Lus5YWI5oqK5a6D5Lus5om+5Ye65p2l77yaCmBgYHtyfQpvdXRsaWVycyA8LSB3aGljaChtazEgIT0gdHJhaW5pbmckRGlyZWN0aW9uKQpsZW5ndGgob3V0bGllcnMpCm91dGxpZXJzWzE6Ml0KYGBgCgrlvILluLjngrnlhbHmnIk0MeS4qu+8jOS4juS4iumdoiBgdGFibGUoKWAg55qE6L6T5Ye6IDIyICsgMTkg5LiA6Ie044CCCuS4i+mdouWIhuaekOS4gOS4i+i/meS6m+W8guW4uOeCueWSjOato+W4uOeCueacieS7gOS5iOS4jeS4gOagt+OAggoK6aaW5YWI5Y+W5YmN5Lik5Liq5byC5bi454K577yI5LiL5qCH5YiG5Yir5Li6MTLlkowxMDjvvInkvZzkuLrliIbmnpDlr7nosaHvvIwK5Y+W5a6D5Lus5ZCO6Z2i55qE5YC85L2c5Li65a+55q+U77yM6KeC5a+f5a6D5Lus5ZCE6Ieq5Zyo5Y2K5b6EMC4wMeiMg+WbtOWGheeahOmCu+Wxhe+8mgpgYGB7cn0KdDEgPC0gMTIKdDIgPC0gMTA4CmVwcyA8LSAwLjAxCmZpbHRlcih0cmFpbmluZywgTGFnMiA+IHRyYWluaW5nJExhZzJbdDFdIC0gZXBzICYgTGFnMiA8IHRyYWluaW5nJExhZzJbdDFdICsgZXBzKQpmaWx0ZXIodHJhaW5pbmcsIExhZzIgPiB0cmFpbmluZyRMYWcyW3QxICsgMV0gLSBlcHMgJiBMYWcyIDwgdHJhaW5pbmckTGFnMlt0MSArIDFdICsgZXBzKQpmaWx0ZXIodHJhaW5pbmcsIExhZzIgPiB0cmFpbmluZyRMYWcyW3QyXSAtIGVwcyAmIExhZzIgPCB0cmFpbmluZyRMYWcyW3QyXSArIGVwcykKZmlsdGVyKHRyYWluaW5nLCBMYWcyID4gdHJhaW5pbmckTGFnMlt0MiArIDFdIC0gZXBzICYgTGFnMiA8IHRyYWluaW5nJExhZzJbdDIgKyAxXSArIGVwcykKYGBgCgrkuI3pmr7nnIvlh7rvvIzov5nkupvngrnkuYvmiYDku6XkuI7igJzoh6rlt7HigJ3nmoQgKkRpcmVjdGlvbiog5Y+W5YC85LiN5ZCM77yM5piv5Zug5Li65pyJ5LiA5Liq5oiW5aSa5Liq55u45ZCM5YC855qE54K55a2Y5Zyo5a+86Ie055qE44CCCgrmnIDlkI7lm57nrZTnrKw46aKY5byV5Ye655qE6Zeu6aKY77yaCgoqIOW9k+aVsOaNrumbhumHjOayoeaciemHjeWPoOeCueaXtu+8jCRLPTEk55qEIGBrbm5gIOWcqOiHqui6q+iuree7g+mbhuS4iueahOato+ehrueOh+aYryAxMDAl77yMCgoqIOWmguaenOaVsOaNrumbhumHjOaciemHjeWPoOeCue+8jOaYr+WQpuWHuuWHuueOsOW8guW4uOeCueWPluWGs+S6jumHjeWPoOeCueeahOagh+iusO+8iOi/memHjOaYryAqRGlyZWN0aW9uKu+8ieWPluWAvOaYr+WQpuS4gOiHtO+8jArlvILluLjngrnotorlpJrvvIzmraPnoa7njofotorkvY7jgIIKCiMjIDEwaAoKTG9naXN0aWMgUmVncmVzc2lvbiBpbiAoZCkgYW5kIExEQSBtZXRob2QgaW4gKGUpIHByb3ZpZGUgdGhlIGJlc3QgcmVzdWx0IChjb3JyZWN0IHJhdGU6IDYyLjUlKSBvbiB0aGlzIGRhdGEuCgojIFF1ZXN0aW9uIDExCgojIyAxMWEKYGBge3J9CmF1dG8xMSA8LSBBdXRvCmF1dG8xMSRtcGcwMSA8LSAwCmF1dG8xMSRtcGcwMVtBdXRvJG1wZyA+IG1lZGlhbihBdXRvJG1wZyldIDwtIDEKYGBgCgojIyAxMWIKYGBge3J9CnBsb3QoYXV0bzExKQpwYXIobWZyb3c9YygzLDMpKQpib3hwbG90KGN5bGluZGVycyB+IG1wZzAxLCBkYXRhID0gYXV0bzExKQpib3hwbG90KGRpc3BsYWNlbWVudCB+IG1wZzAxLCBkYXRhID0gYXV0bzExKQpib3hwbG90KGhvcnNlcG93ZXIgfiBtcGcwMSwgZGF0YSA9IGF1dG8xMSkKYm94cGxvdCh3ZWlnaHQgfiBtcGcwMSwgZGF0YSA9IGF1dG8xMSkKYm94cGxvdChhY2NlbGVyYXRpb24gfiBtcGcwMSwgZGF0YSA9IGF1dG8xMSkKYm94cGxvdCh5ZWFyIH4gbXBnMDEsIGRhdGEgPSBhdXRvMTEpCmJveHBsb3Qob3JpZ2luIH4gbXBnMDEsIGRhdGEgPSBhdXRvMTEpCnBhcihtZnJvdyA9IGMoMSwxKSkKY29yKGF1dG8xMVssIC05XSkgICMgY2FsY3VsYXRlIGNvciBtYXRyaXggZXhjZXB0ICpuYW1lKiBjb2x1bW4KYGBgCgoqY3lsaW5kZXIqLCAqZGlzcGxhY2VtZW50KiwgKmhvcnNlcG93ZXIqIGFuZCAqd2VpZ2h0KiBoYXZlIHJlbGF0aXZlIGhpZ2ggbmVnYXRpdmUgY29yZWxhdGlvbnMgd2l0aCAqbXBnMDEqIChhYnNvbHV0ZSB2YWx1ZSBhcmUgYWxsIGFib3ZlIDAuNjUpLgoKIyMgMTFjCmBgYHtyfQppc1RyYWluaW5nIDwtIHNhbXBsZShucm93KGF1dG8xMSksIG5yb3coYXV0bzExKSAqIDAuOCkKdHJhaW5pbmcgPC0gYXV0bzExW2lzVHJhaW5pbmcsIF0KdGVzdGluZyA8LSBhdXRvMTFbLWlzVHJhaW5pbmcsIF0KYGBgCgojIyAxMWQKYGBge3J9CmYxMWQgPC0gbGRhKG1wZzAxIH4gY3lsaW5kZXJzICsgZGlzcGxhY2VtZW50ICsgaG9yc2Vwb3dlciArIHdlaWdodCwgZGF0YSA9IGF1dG8xMSkKcHJlZDExZCA8LSBwcmVkaWN0KGYxMWQsIHRlc3RpbmcpCm1lYW4ocHJlZDExZCRjbGFzcyA9PSB0ZXN0aW5nJG1wZzAxKQpgYGAKVGhlIHRlc3QgZXJyb3IgaXMgJDEgLSAwLjg3MyA9IDAuMTI3JSQuCgojIyAxMWUKYGBge3J9CmYxMWUgPC0gcWRhKG1wZzAxIH4gY3lsaW5kZXJzICsgZGlzcGxhY2VtZW50ICsgaG9yc2Vwb3dlciArIHdlaWdodCwgZGF0YSA9IGF1dG8xMSkKcHJlZDExZSA8LSBwcmVkaWN0KGYxMWUsIHRlc3RpbmcpCjEgLSBtZWFuKHByZWQxMWUkY2xhc3MgPT0gdGVzdGluZyRtcGcwMSkKYGBgCgojIyAxMWYKYGBge3J9CmYxMWYgPC0gZ2xtKG1wZzAxIH4gY3lsaW5kZXJzICsgZGlzcGxhY2VtZW50ICsgaG9yc2Vwb3dlciArIHdlaWdodCwgZGF0YSA9IGF1dG8xMSkKcHJlZDExZiA8LSBwcmVkaWN0KGYxMWYsIHRlc3RpbmcpCnJlczExZiA8LSByZXAoMCwgbGVuZ3RoKHByZWQxMWYpKQpyZXMxMWZbcHJlZDExZiA+IDAuNV0gPSAxCjEgLSBtZWFuKHJlczExZiA9PSB0ZXN0aW5nJG1wZzAxKQpgYGAKCiMjIDExZwpgYGB7cn0Kc2V0LnNlZWQoMSkKdHJhaW4xMWcgPC0gYXV0bzExW2lzVHJhaW5pbmcsIGMoJ2N5bGluZGVycycsICdkaXNwbGFjZW1lbnQnLCAnaG9yc2Vwb3dlcicsICd3ZWlnaHQnKV0KdGVzdDExZyA8LSBhdXRvMTFbLWlzVHJhaW5pbmcsIGMoJ2N5bGluZGVycycsICdkaXNwbGFjZW1lbnQnLCAnaG9yc2Vwb3dlcicsICd3ZWlnaHQnKV0KdHJhaW5SZXMgPC0gYXV0bzExJG1wZzAxW2lzVHJhaW5pbmddCmYxMWcxIDwtIGtubih0cmFpbjExZywgdGVzdDExZywgdHJhaW5SZXMsIGsgPSAxKQptZWFuKGYxMWcxID09IGF1dG8xMSRtcGcwMVstaXNUcmFpbmluZ10pCmYxMWczIDwtIGtubih0cmFpbjExZywgdGVzdDExZywgdHJhaW5SZXMsIGsgPSAzKQptZWFuKGYxMWczID09IGF1dG8xMSRtcGcwMVstaXNUcmFpbmluZ10pCmYxMWc1IDwtIGtubih0cmFpbjExZywgdGVzdDExZywgdHJhaW5SZXMsIGsgPSA1KQptZWFuKGYxMWc1ID09IGF1dG8xMSRtcGcwMVstaXNUcmFpbmluZ10pCmYxMWcxMDAgPC0ga25uKHRyYWluMTFnLCB0ZXN0MTFnLCB0cmFpblJlcywgayA9IDEwMCkKbWVhbihmMTFnMTAwID09IGF1dG8xMSRtcGcwMVstaXNUcmFpbmluZ10pCmBgYApUaGUgbWluaW11bSBlcnJvciByYXRlIGlzICQxIC0gMC45MzcgPSAwLjA2MyQuCiRrID0gNSQgcGVyZm9ybXMgYmVzdCBvbiB0aGlzIGRhdGEgc2V0LgoKIyBRdWVzdGlvbiAxMgojIyAxMmEKYGBge3J9ClBvd2VyIDwtIGZ1bmN0aW9uKCkgeyAyIF4gMyB9ClBvd2VyKCkKYGBgCgojIyAxMmIgJiAxMmMKYGBge3J9ClBvd2VyMiA8LSBmdW5jdGlvbih4LCBhKSB7IHggXiBhIH0KUG93ZXIyKDMsIDgpClBvd2VyMigxMCwgMykKUG93ZXIyKDgsIDE3KQpQb3dlcjIoMTMxLCAzKQpgYGAKCiMjIDEyZApBcyB0aGUgY3VycmVudCB2ZXJzaW9uIG9mIFIsIHRoZSBsYXN0IGV4cHJlc3Npb24gaXMgdGhlIHJldHVybiB2YWx1ZSBvZiB0aGUgZnVuY3Rpb24uClNvIGl0J3MgdW5uZWNlc3NhcnkgdG8gdXNlIGByZXR1cm4oKWAgZnVuY3Rpb24gZXhwbGljaXRseToKYGBge3J9CnZlcnNpb24KYGBgCgojIyAxMmUKYGBge3J9CnBhcihtZnJvdz1jKDEsMikpCnggPC0gMToxMApwbG90KHgsIFBvd2VyMih4LCAyKSwgbWFpbiA9ICdmKHgpID0geF4yJywgeGxhYiA9ICd4ID0gMToxMCcsIHlsYWIgPSAneF4yJykKcGxvdChsb2coeCksIGxvZyhQb3dlcjIoeCwgMikpLCBtYWluID0gJ2xvZyhmKHgpKScsIHhsYWIgPSAneCA9IDE6MTAnLCB5bGFiID0gJzIgKiBsb2coeCknKQpgYGAKCiMjIDEyZgpgYGB7cn0KUGxvdFBvd2VyIDwtIGZ1bmN0aW9uKHgsIHApIHsKICBwbG90KHgsIHggXiBwLCBtYWluID0gJ3Bvd2VyIHAgb2YgeCcpCn0KUGxvdFBvd2VyKDE6MTAsIDMpCmBgYAoKIyBRdWVzdGlvbiAxMwpFeHRyYWN0IHRyYWluaW5nIGFuZCB0ZXN0aW5nIGRhdGEsIGFuZCBwcmVkaWN0IGJhc2VkIG9uIHRoZSBmaXJzdCBncm91cCBvZiBwcmVkaWN0b3JzIHdpdGggbG9naXN0aWMgcmVncmVzc2lvbiBtZXRob2Q6CmBgYHtyfQpib3N0b24xMyA8LSBCb3N0b24KYm9zdG9uMTMkY3JpbWUucmF0ZSA8LSAwCmJvc3RvbjEzJGNyaW1lLnJhdGVbQm9zdG9uJGNyaW0gPiBtZWRpYW4oQm9zdG9uJGNyaW0pXSA8LSAxCmlzVHJhaW5pbmcgPC0gc2FtcGxlKG5yb3coYm9zdG9uMTMpLCAwLjggKiBucm93KGJvc3RvbjEzKSkKdHJhaW5pbmcgPC0gYm9zdG9uMTNbaXNUcmFpbmluZywgXQp0ZXN0aW5nIDwtIGJvc3RvbjEzWy1pc1RyYWluaW5nLCBdCmxyLmZpdDEgPC0gZ2xtKGNyaW1lLnJhdGUgfiB6biArIGluZHVzICsgY2hhcyArIG5veCwgZGF0YSA9IGJvc3RvbjEzLCBmYW1pbHkgPSAnYmlub21pYWwnKQpsci5wcmVkMSA8LSBwcmVkaWN0KGxyLmZpdDEsIHRlc3RpbmcpCmxyLnJlczEgPC0gcmVwKDAsIGxlbmd0aChsci5wcmVkMSkpCmxyLnJlczFbbHIucHJlZDEgPiBtZWRpYW4oQm9zdG9uJGNyaW0pXSA8LSAxCm1lYW4obHIucmVzMSA9PSB0ZXN0aW5nJGNyaW1lLnJhdGUpCmBgYAoKVXNpbmcgYW5vdGhlciBncm91cCBvZiBwcmVkaWN0b3JzOgpgYGB7cn0KbHIuZml0MiA8LSBnbG0oY3JpbWUucmF0ZSB+IHJtICsgYWdlICsgZGlzICsgcmFkLCBkYXRhID0gYm9zdG9uMTMsIGZhbWlseSA9ICdiaW5vbWlhbCcpCmxyLnByZWQyIDwtIHByZWRpY3QobHIuZml0MiwgdGVzdGluZykKbHIucmVzMiA8LSByZXAoMCwgbGVuZ3RoKGxyLnByZWQyKSkKbHIucmVzMltsci5wcmVkMiA+IG1lZGlhbihCb3N0b24kY3JpbSldIDwtIDEKbWVhbihsci5yZXMyID09IHRlc3RpbmckY3JpbWUucmF0ZSkKYGBgCgpQcmVkaWN0IG9uIHRoZSAxc3QgZ3JvdXAgcHJlZGljdG9ycyB3aXRoIExEQToKYGBge3J9CmxkYS5maXQxIDwtIGxkYShjcmltZS5yYXRlIH4gem4gKyBpbmR1cyArIGNoYXMgKyBub3gsIGRhdGEgPSBib3N0b24xMykKbGRhLnByZWQxIDwtIHByZWRpY3QobGRhLmZpdDEsIHRlc3RpbmcpCm1lYW4obGRhLnByZWQxJGNsYXNzID09IHRlc3RpbmckY3JpbWUucmF0ZSkKYGBgCgpQcmVkaWN0IG9uIHRoZSAybmQgZ3JvdXAgcHJlZGljdG9ycyB3aXRoIExEQToKYGBge3J9CmxkYS5maXQyIDwtIGxkYShjcmltZS5yYXRlIH4gcm0gKyBhZ2UgKyBkaXMgKyByYWQsIGRhdGEgPSBib3N0b24xMykKbGRhLnByZWQyIDwtIHByZWRpY3QobGRhLmZpdDIsIHRlc3RpbmcpCm1lYW4obGRhLnByZWQyJGNsYXNzID09IHRlc3RpbmckY3JpbWUucmF0ZSkKYGBgCgpQcmVkaWN0IG9uIHRoZSAxc3QgZ3JvdXAgcHJlZGljdG9ycyB3aXRoIEtOTiAoayA9IDUpOgpgYGB7cn0Ka3RyYWluMSA8LSBib3N0b24xM1tpc1RyYWluaW5nLCBjKCd6bicsICdpbmR1cycsICdjaGFzJywgJ25veCcpXQprdGVzdDEgPC0gYm9zdG9uMTNbLWlzVHJhaW5pbmcsIGMoJ3puJywgJ2luZHVzJywgJ2NoYXMnLCAnbm94JyldCmt0cmFpbi5yZXN1bHQgPC0gYm9zdG9uMTMkY3JpbWUucmF0ZVtpc1RyYWluaW5nXQprbW9kMSA8LSBrbm4oa3RyYWluMSwga3Rlc3QxLCBrdHJhaW4ucmVzdWx0LCBrID0gNSkKbWVhbihrbW9kMSA9PSBib3N0b24xMyRjcmltZS5yYXRlWy1pc1RyYWluaW5nXSkKYGBgCgpQcmVkaWN0IG9uIHRoZSAybmQgZ3JvdXAgcHJlZGljdG9ycyB3aXRoIEtOTiAoayA9IDUpOgpgYGB7cn0Ka3RyYWluMiA8LSBib3N0b24xM1tpc1RyYWluaW5nLCBjKCdybScsICdhZ2UnLCAnZGlzJywgJ3JhZCcpXQprdGVzdDIgPC0gYm9zdG9uMTNbLWlzVHJhaW5pbmcsIGMoJ3JtJywgJ2FnZScsICdkaXMnLCAncmFkJyldCmttb2QyIDwtIGtubihrdHJhaW4yLCBrdGVzdDIsIGt0cmFpbi5yZXN1bHQsIGsgPSA1KQptZWFuKGttb2QyID09IGJvc3RvbjEzJGNyaW1lLnJhdGVbLWlzVHJhaW5pbmddKQpgYGAKCktOTiAoaz01KSBiYXNlZCBvbiB0aGUgZmlyc3QgZ3JvdXAgb2YgcHJlZGljdG9ycyBoYXMgdGhlIGhpZ2hlc3QgYWNjdXJhY3ku