5.3.1 The validation Set Approach
Se usa para estimar las tasas de error que resultan despues de fijar varios modelos lineales en el dataset. Se utilizará el dataset Auto
.
Antes de comenzar se usara la funcion set.seed()
para fijar el numero random del generador, esto es para que los resultados concuerden con los fijados en el laboratorio guiado del libro ISLR. Se utilizara la funcion sample()
para dividir el dataset en dos partes, tomando estas observaciones como el set de entrenamiento.
library(ISLR)
set.seed(1)
train <- sample(392, 196)
Se usa el parametro subset
en la funcion lm()
para usar solamente las observaciones correspondientes al training set.
lm.fit <- lm(mpg ~ horsepower, data = Auto, subset = train)
Ahora se usa la funcion predict()
para estimar la respuesta de todas las 392 observaciones y obtenemos el MSE de las 196 observaciones en el set de validacion.
attach(Auto)
mean((mpg - predict(lm.fit, Auto))[-train]^2)
[1] 26.14142
Por lo tanto el MSE del test estimado para la regresion lineal es \(26.14\). Ahora usamos la funcion poly()
para estimar el error del test para las regresiones polinomiales.
lm.fit2 <- lm(mpg ~ poly(horsepower,2), data = Auto, subset = train)
mean((mpg - predict(lm.fit2, Auto))[-train]^2)
[1] 19.82259
lm.fit3 <- lm(mpg ~ poly(horsepower,3), data = Auto, subset = train)
mean((mpg - predict(lm.fit3, Auto))[-train]^2)
[1] 19.78252
Las tasas de eroor son \(19.82\) y \(19.78\) respectivamente. Si se usa un training set distinto, entonces se van a obtener errores diferentes en el set de validacion.
set.seed(2)
train <- sample(392, 196)
lm.fit <- lm(mpg ~ horsepower, subset = train)
mean((mpg - predict(lm.fit, Auto))[-train]^2)
[1] 23.29559
lm.fit2 <- lm(mpg ~ poly(horsepower,2), data = Auto, subset = train)
mean((mpg - predict(lm.fit2, Auto))[-train]^2)
[1] 18.90124
lm.fit3 <- lm(mpg ~ poly(horsepower,3), data = Auto, subset = train)
mean((mpg - predict(lm.fit3, Auto))[-train]^2)
[1] 19.2574
Estos resultados son consistentes con los descubrimientos anteriores. Un modelo que predice mpg
usando una funcion cuadratica de los horsepower
se desempeña mejor que un modelo que involucra solo una funcion lineal de horsepower
, y hay una minima evidencia a favor del modelo que usa una funcion cubica en horsepower
.
Leave-One-Out Cross-Validation
El estimador LOOCV puede ser computado automaticamente para cualquier modelo lineal generalizado usando glm()
y cv.glm()
. La funcion glm()
se uso para desarrollar una regresion logistica enviandole como argumento family=binomial
, pero si se utiliza sin ese argumento, se desarrollara una regresion lineal como la lm()
:
glm.fit <- glm(mpg ~ horsepower, data = Auto)
coef(glm.fit)
(Intercept) horsepower
39.9358610 -0.1578447
lm.fit <- lm(mpg ~ horsepower, data = Auto)
coef(lm.fit)
(Intercept) horsepower
39.9358610 -0.1578447
Por lo que se puede ver que son modelos de regresion lineal identicos. Para el laboratorio se utilizara la funcion glm()
para obtener la regresion lineal, ya que la funcion cv.glm()
puede ser utilizada, esta funcion es parte de la libreria boot
.
library(boot)
glm.fit <- glm(mpg ~ horsepower, data = Auto)
cv.err <- cv.glm(Auto, glm.fit)
cv.err$delta
[1] 24.23151 24.23114
La funcion cv.glm()
produce una lista con distintos componentes. Los dos numeros en el vector delta contienen los resultados de la cross-validation. En este caso los numeros son identicos, y corresponen al estadistico LOOCV.
Ahora se va a repetir este procedimiento para nuevos modelos incrementando la complejidad polinomial. Se usara la funcion for()
para automatizar el proceso, los errores se iran guardando en el vector cv.error
.
cv.error <- rep(0,5)
for(i in 1:5){
glm.fit <- glm(mpg ~ poly(horsepower, i), data = Auto)
cv.error[i] <- cv.glm(Auto, glm.fit)$delta[1]
}
cv.error
[1] 24.23151 19.24821 19.33498 19.42443 19.03321
como vimos anteriormente el MSE estimado del test entre los modelos lineales y cuadraticos hay una disminucion considerable, pero no hay una mejora al usar polinomios mas altos.
k-Fold Cross-Validation
La funcion cv.glm()
tambien puede ser usanda para implementar el k-fold CV. Usaremos 10 folds en el dataset Auto.
set.seed(17)
cv.error.10 <- rep(0,10)
for(i in 1:10){
glm.fit <- glm(mpg ~ poly(horsepower,i), data = Auto)
cv.error.10[i] <- cv.glm(Auto,glm.fit, K = 10)$delta[1]
}
cv.error.10
[1] 24.20520 19.18924 19.30662 19.33799 18.87911 19.02103
[7] 18.89609 19.71201 18.95140 19.50196
The Bootstrap
Estimating the Accuracy of a Statistic of Interest
Una de las grandes ventajas del enfoque bootstrap es que puede ser aplicado en casi todas las situaciones. Los calculos matematicos complicados no son requeridos. Realizar un analisis bootstrap en R conlleva solo dos pasos.
- Primero se debe crear una funcion que compute el estadistico de interes.
- Segundo usamos la funcion
boot()
(de la libreria boot
) para realizar el bootstrap.
El dataset Portafolio
se utilizara en este ejemplo. Para ilustrar el uso de bootstrap en este dataset primero se debe crear una funcion alpha.fn()
que toma \((X,Y)\) como input como un vector indicando cuales observaciones deberian ser usadas para estimar \(\alpha\).
alpha.fn <- function(data,index){
X <- data$X[index]
Y <- data$Y[index]
return((var(Y) - cov(X,Y)) / (var(X) + var(Y) - 2*cov(X,Y)))
}
Esta funcion devuelve como output un estimado para \(\alpha\) basado en el argumento index
. Ahora se usara para estimar \(\alpha\) usando 100 observaciones:
alpha.fn(Portfolio, 1:100)
[1] 0.5758321
Ahora se hara un sample()
para seleccionar 100 observaciones del rango 1 a 100.
set.seed(1)
alpha.fn(Portfolio, sample(100,100,replace = TRUE))
[1] 0.5963833
Se puede implementar un analisis bootstrap realizando este comando muchas veces. Esto lo puede hacer la funcion boot()
, se produciran \(R=1000\) estimaciones bootstrap para \(\alpha\).
boot(Portfolio, alpha.fn, R = 1000)
ORDINARY NONPARAMETRIC BOOTSTRAP
Call:
boot(data = Portfolio, statistic = alpha.fn, R = 1000)
Bootstrap Statistics :
original bias std. error
t1* 0.5758321 -7.315422e-05 0.08861826
El output muestra que usando el dataset original \(\hat{\alpha} = 0.5758\) y que el bootstrap estimado es \(SE(\hat{\alpha}= 0.0886\).
Estimating the Accuracy of a Linear Regression Model
El enfoque bootstrap puede ser usado para evaluar la variablidad de las estimaciones de los coeficientes y las predicciones del metodo de aprendizaje estadistico. Se utilizara bootstrap para evaluar la variabilidad de los estimados para \(\beta_0\) y \(\beta_1\), el intercepto y la pendiente para el modelo de regresion lineal que usa horsepower
para predecir mpg
en el dataset Auto
.
Primero se va a crear una funcion simpre boot.fn()
, la cual toma el dataset Auto
y los indices de observacion y regresa el intercepto y la pendiente estimadas para el modelo de regresion lineal. Despues aplicamos esta funcion al dataset completo y computaremos los esmimados de \(\beta_0\) y de \(\beta_1\).
boot.fn <- function(data, index)
return(coef(lm(mpg ~ horsepower, data = data, subset = index)))
boot.fn(Auto, 1:392)
(Intercept) horsepower
39.9358610 -0.1578447
La funcion boot.fn()
puede ser usara para crear estimadores bootstrap para los interceptos y la pendiente haciendo un muestreo de las observaciones:
set.seed(1)
boot.fn(Auto, sample(392, 392, replace = TRUE))
(Intercept) horsepower
38.7387134 -0.1481952
boot.fn(Auto,sample(392, 392, replace = TRUE))
(Intercept) horsepower
40.0383086 -0.1596104
Ahora usaremos la funcion boot()
para obtener los errores estandar de \(1000\) estimaciones bootstrap para el intercepto y la pendiente.
boot(Auto, boot.fn, 1000)
ORDINARY NONPARAMETRIC BOOTSTRAP
Call:
boot(data = Auto, statistic = boot.fn, R = 1000)
Bootstrap Statistics :
original bias std. error
t1* 39.9358610 0.02972191 0.860007896
t2* -0.1578447 -0.00030823 0.007404467
Esto indica que el estimador bootstrap es \(SE(\hat{\beta_0}) = 0.86\) y \(SE(\hat{\beta_1}) = 0.0074\).
Las formulas estardar pueden ser usadas para obtener los erroes estandar de los coeficientes regresores en un modelo lineal:
summary(lm(mpg ~ horsepower, data = Auto))$coef
Estimate Std. Error t value Pr(>|t|)
(Intercept) 39.9358610 0.717498656 55.65984 1.220362e-187
horsepower -0.1578447 0.006445501 -24.48914 7.031989e-81
boot.fn <- function(data, index)
coefficients(lm(mpg ~ horsepower + I(horsepower^2), data = data, subset = index))
set.seed(1)
boot(Auto, boot.fn, 1000)
ORDINARY NONPARAMETRIC BOOTSTRAP
Call:
boot(data = Auto, statistic = boot.fn, R = 1000)
Bootstrap Statistics :
original bias std. error
t1* 56.900099702 6.098115e-03 2.0944855842
t2* -0.466189630 -1.777108e-04 0.0334123802
t3* 0.001230536 1.324315e-06 0.0001208339
summary(lm( mpg ~ horsepower + I(horsepower^2), data = Auto))$coef
Estimate Std. Error t value Pr(>|t|)
(Intercept) 56.900099702 1.8004268063 31.60367 1.740911e-109
horsepower -0.466189630 0.0311246171 -14.97816 2.289429e-40
I(horsepower^2) 0.001230536 0.0001220759 10.08009 2.196340e-21
LS0tDQp0aXRsZTogIjUuMyBMYWI6IENyb3NzLVZhbGlkYXRpb24gYW5kIHRoZSBCb290c3RyYXAiDQpvdXRwdXQ6IGh0bWxfbm90ZWJvb2sNCi0tLQ0KDQojIyBQYXVsYSBDYXphbGkNCiMjIyBGaWFiaWxpZGFkDQoNCiMjIDUuMy4xIFRoZSB2YWxpZGF0aW9uIFNldCBBcHByb2FjaA0KDQpTZSB1c2EgcGFyYSBlc3RpbWFyIGxhcyB0YXNhcyBkZSBlcnJvciBxdWUgcmVzdWx0YW4gZGVzcHVlcyBkZSBmaWphciB2YXJpb3MgbW9kZWxvcyBsaW5lYWxlcyBlbiBlbCBkYXRhc2V0LiANClNlIHV0aWxpemFyw6EgZWwgZGF0YXNldCBgQXV0b2AuDQoNCkFudGVzIGRlIGNvbWVuemFyIHNlIHVzYXJhIGxhIGZ1bmNpb24gYHNldC5zZWVkKClgIHBhcmEgZmlqYXIgZWwgbnVtZXJvIHJhbmRvbSBkZWwgZ2VuZXJhZG9yLCBlc3RvIGVzIHBhcmEgcXVlIGxvcyByZXN1bHRhZG9zIGNvbmN1ZXJkZW4gY29uIGxvcyBmaWphZG9zIGVuIGVsIGxhYm9yYXRvcmlvIGd1aWFkbyBkZWwgbGlicm8gSVNMUi4gDQpTZSB1dGlsaXphcmEgbGEgZnVuY2lvbiBgc2FtcGxlKClgIHBhcmEgZGl2aWRpciBlbCBkYXRhc2V0IGVuIGRvcyBwYXJ0ZXMsIHRvbWFuZG8gZXN0YXMgb2JzZXJ2YWNpb25lcyBjb21vIGVsIHNldCBkZSBlbnRyZW5hbWllbnRvLg0KDQpgYGB7cn0NCmxpYnJhcnkoSVNMUikNCnNldC5zZWVkKDEpDQp0cmFpbiA8LSBzYW1wbGUoMzkyLCAxOTYpDQpgYGANCg0KU2UgdXNhIGVsIHBhcmFtZXRybyBgc3Vic2V0YCBlbiBsYSBmdW5jaW9uIGBsbSgpYCBwYXJhIHVzYXIgc29sYW1lbnRlIGxhcyBvYnNlcnZhY2lvbmVzIGNvcnJlc3BvbmRpZW50ZXMgYWwgdHJhaW5pbmcgc2V0Lg0KYGBge3J9DQpsbS5maXQgPC0gbG0obXBnIH4gaG9yc2Vwb3dlciwgZGF0YSA9IEF1dG8sIHN1YnNldCA9IHRyYWluKQ0KYGBgDQoNCkFob3JhIHNlIHVzYSBsYSBmdW5jaW9uIGBwcmVkaWN0KClgIHBhcmEgZXN0aW1hciBsYSByZXNwdWVzdGEgZGUgdG9kYXMgbGFzIDM5MiBvYnNlcnZhY2lvbmVzIHkgb2J0ZW5lbW9zIGVsIE1TRSBkZSBsYXMgMTk2IG9ic2VydmFjaW9uZXMgZW4gZWwgc2V0IGRlIHZhbGlkYWNpb24uIA0KYGBge3J9DQphdHRhY2goQXV0bykNCmBgYA0KYGBge3J9DQptZWFuKChtcGcgLSBwcmVkaWN0KGxtLmZpdCwgQXV0bykpWy10cmFpbl1eMikNCmBgYA0KDQpQb3IgbG8gdGFudG8gZWwgTVNFIGRlbCB0ZXN0IGVzdGltYWRvIHBhcmEgbGEgcmVncmVzaW9uIGxpbmVhbCBlcyAkMjYuMTQkLg0KQWhvcmEgdXNhbW9zIGxhIGZ1bmNpb24gYHBvbHkoKWAgcGFyYSBlc3RpbWFyIGVsIGVycm9yIGRlbCB0ZXN0IHBhcmEgbGFzIHJlZ3Jlc2lvbmVzIHBvbGlub21pYWxlcy4NCmBgYHtyfQ0KbG0uZml0MiA8LSBsbShtcGcgfiBwb2x5KGhvcnNlcG93ZXIsMiksIGRhdGEgPSBBdXRvLCBzdWJzZXQgPSB0cmFpbikNCm1lYW4oKG1wZyAtIHByZWRpY3QobG0uZml0MiwgQXV0bykpWy10cmFpbl1eMikNCmxtLmZpdDMgPC0gbG0obXBnIH4gcG9seShob3JzZXBvd2VyLDMpLCBkYXRhID0gQXV0bywgc3Vic2V0ID0gdHJhaW4pDQptZWFuKChtcGcgLSBwcmVkaWN0KGxtLmZpdDMsIEF1dG8pKVstdHJhaW5dXjIpDQpgYGANCg0KTGFzIHRhc2FzIGRlIGVyb29yIHNvbiAkMTkuODIkIHkgJDE5Ljc4JCByZXNwZWN0aXZhbWVudGUuIFNpIHNlIHVzYSB1biB0cmFpbmluZyBzZXQgZGlzdGludG8sIGVudG9uY2VzIHNlIHZhbiBhIG9idGVuZXIgZXJyb3JlcyBkaWZlcmVudGVzIGVuIGVsIHNldCBkZSB2YWxpZGFjaW9uLg0KYGBge3J9DQpzZXQuc2VlZCgyKQ0KdHJhaW4gPC0gc2FtcGxlKDM5MiwgMTk2KQ0KbG0uZml0IDwtIGxtKG1wZyB+IGhvcnNlcG93ZXIsIHN1YnNldCA9IHRyYWluKQ0KbWVhbigobXBnIC0gcHJlZGljdChsbS5maXQsIEF1dG8pKVstdHJhaW5dXjIpDQpsbS5maXQyIDwtIGxtKG1wZyB+IHBvbHkoaG9yc2Vwb3dlciwyKSwgZGF0YSA9IEF1dG8sIHN1YnNldCA9IHRyYWluKQ0KbWVhbigobXBnIC0gcHJlZGljdChsbS5maXQyLCBBdXRvKSlbLXRyYWluXV4yKQ0KbG0uZml0MyA8LSBsbShtcGcgfiBwb2x5KGhvcnNlcG93ZXIsMyksIGRhdGEgPSBBdXRvLCBzdWJzZXQgPSB0cmFpbikNCm1lYW4oKG1wZyAtIHByZWRpY3QobG0uZml0MywgQXV0bykpWy10cmFpbl1eMikNCg0KYGBgDQoNCkVzdG9zIHJlc3VsdGFkb3Mgc29uIGNvbnNpc3RlbnRlcyBjb24gbG9zIGRlc2N1YnJpbWllbnRvcyBhbnRlcmlvcmVzLiBVbiBtb2RlbG8gcXVlIHByZWRpY2UgYG1wZ2AgdXNhbmRvIHVuYSBmdW5jaW9uIGN1YWRyYXRpY2EgZGUgbG9zIGBob3JzZXBvd2VyYCBzZSBkZXNlbXBlw7FhIG1lam9yIHF1ZSB1biBtb2RlbG8gcXVlIGludm9sdWNyYSBzb2xvIHVuYSBmdW5jaW9uIGxpbmVhbCBkZSBgaG9yc2Vwb3dlcmAsIHkgaGF5IHVuYSBtaW5pbWEgZXZpZGVuY2lhIGEgZmF2b3IgZGVsIG1vZGVsbyBxdWUgdXNhIHVuYSBmdW5jaW9uIGN1YmljYSBlbiBgaG9yc2Vwb3dlcmAuDQoNCiMjIExlYXZlLU9uZS1PdXQgQ3Jvc3MtVmFsaWRhdGlvbg0KDQpFbCBlc3RpbWFkb3IgTE9PQ1YgIHB1ZWRlIHNlciBjb21wdXRhZG8gYXV0b21hdGljYW1lbnRlIHBhcmEgY3VhbHF1aWVyIG1vZGVsbyBsaW5lYWwgZ2VuZXJhbGl6YWRvIHVzYW5kbyBgZ2xtKClgIHkgYGN2LmdsbSgpYC4gTGEgZnVuY2lvbiBgZ2xtKClgIHNlIHVzbyBwYXJhIGRlc2Fycm9sbGFyIHVuYSByZWdyZXNpb24gbG9naXN0aWNhIGVudmlhbmRvbGUgY29tbyBhcmd1bWVudG8gYGZhbWlseT1iaW5vbWlhbGAsIHBlcm8gc2kgc2UgdXRpbGl6YSBzaW4gZXNlIGFyZ3VtZW50bywgc2UgZGVzYXJyb2xsYXJhIHVuYSByZWdyZXNpb24gbGluZWFsIGNvbW8gbGEgYGxtKClgOg0KYGBge3J9DQpnbG0uZml0IDwtIGdsbShtcGcgfiBob3JzZXBvd2VyLCBkYXRhID0gQXV0bykNCmNvZWYoZ2xtLmZpdCkNCmBgYA0KDQpgYGB7cn0NCmxtLmZpdCA8LSBsbShtcGcgfiBob3JzZXBvd2VyLCBkYXRhID0gQXV0bykNCmNvZWYobG0uZml0KQ0KYGBgDQoNClBvciBsbyBxdWUgc2UgcHVlZGUgdmVyIHF1ZSBzb24gbW9kZWxvcyBkZSByZWdyZXNpb24gbGluZWFsIGlkZW50aWNvcy4gUGFyYSBlbCBsYWJvcmF0b3JpbyBzZSB1dGlsaXphcmEgbGEgZnVuY2lvbiBgZ2xtKClgIHBhcmEgb2J0ZW5lciBsYSByZWdyZXNpb24gbGluZWFsLCB5YSBxdWUgbGEgZnVuY2lvbiBgY3YuZ2xtKClgIHB1ZWRlIHNlciB1dGlsaXphZGEsIGVzdGEgZnVuY2lvbiBlcyBwYXJ0ZSBkZSBsYSBsaWJyZXJpYSBgYm9vdGAuDQpgYGB7cn0NCmxpYnJhcnkoYm9vdCkNCmdsbS5maXQgPC0gZ2xtKG1wZyB+IGhvcnNlcG93ZXIsIGRhdGEgPSBBdXRvKQ0KY3YuZXJyIDwtIGN2LmdsbShBdXRvLCBnbG0uZml0KQ0KY3YuZXJyJGRlbHRhDQpgYGANCg0KTGEgZnVuY2lvbiBgY3YuZ2xtKClgIHByb2R1Y2UgdW5hIGxpc3RhIGNvbiBkaXN0aW50b3MgY29tcG9uZW50ZXMuIExvcyBkb3MgbnVtZXJvcyBlbiBlbCB2ZWN0b3IgKmRlbHRhKiBjb250aWVuZW4gbG9zIHJlc3VsdGFkb3MgZGUgbGEgY3Jvc3MtdmFsaWRhdGlvbi4gRW4gZXN0ZSBjYXNvIGxvcyBudW1lcm9zIHNvbiBpZGVudGljb3MsIHkgY29ycmVzcG9uZW4gYWwgZXN0YWRpc3RpY28gTE9PQ1YuIA0KDQpBaG9yYSBzZSB2YSBhIHJlcGV0aXIgZXN0ZSBwcm9jZWRpbWllbnRvIHBhcmEgbnVldm9zIG1vZGVsb3MgaW5jcmVtZW50YW5kbyBsYSBjb21wbGVqaWRhZCBwb2xpbm9taWFsLiBTZSB1c2FyYSBsYSBmdW5jaW9uIGBmb3IoKWAgcGFyYSBhdXRvbWF0aXphciBlbCBwcm9jZXNvLCBsb3MgZXJyb3JlcyBzZSBpcmFuIGd1YXJkYW5kbyBlbiBlbCB2ZWN0b3IgYGN2LmVycm9yYC4NCmBgYHtyfQ0KY3YuZXJyb3IgPC0gcmVwKDAsNSkNCmZvcihpIGluIDE6NSl7DQogIGdsbS5maXQgPC0gZ2xtKG1wZyB+IHBvbHkoaG9yc2Vwb3dlciwgaSksIGRhdGEgPSBBdXRvKQ0KICBjdi5lcnJvcltpXSA8LSBjdi5nbG0oQXV0bywgZ2xtLmZpdCkkZGVsdGFbMV0NCn0NCmN2LmVycm9yDQpgYGANCg0KY29tbyB2aW1vcyBhbnRlcmlvcm1lbnRlIGVsIE1TRSBlc3RpbWFkbyBkZWwgdGVzdCBlbnRyZSBsb3MgbW9kZWxvcyBsaW5lYWxlcyB5IGN1YWRyYXRpY29zIGhheSB1bmEgZGlzbWludWNpb24gY29uc2lkZXJhYmxlLCBwZXJvIG5vIGhheSB1bmEgbWVqb3JhIGFsIHVzYXIgcG9saW5vbWlvcyBtYXMgYWx0b3MuDQoNCiMjIGstRm9sZCBDcm9zcy1WYWxpZGF0aW9uDQoNCkxhIGZ1bmNpb24gYGN2LmdsbSgpYCB0YW1iaWVuIHB1ZWRlIHNlciB1c2FuZGEgcGFyYSBpbXBsZW1lbnRhciBlbCBrLWZvbGQgQ1YuIFVzYXJlbW9zIDEwIGZvbGRzIGVuIGVsIGRhdGFzZXQgQXV0by4gDQpgYGB7cn0NCnNldC5zZWVkKDE3KQ0KY3YuZXJyb3IuMTAgPC0gcmVwKDAsMTApDQpmb3IoaSBpbiAxOjEwKXsNCiAgZ2xtLmZpdCA8LSBnbG0obXBnIH4gcG9seShob3JzZXBvd2VyLGkpLCBkYXRhID0gQXV0bykNCiAgY3YuZXJyb3IuMTBbaV0gPC0gY3YuZ2xtKEF1dG8sZ2xtLmZpdCwgSyA9IDEwKSRkZWx0YVsxXQ0KfQ0KY3YuZXJyb3IuMTANCmBgYA0KDQojIyBUaGUgQm9vdHN0cmFwDQojIyMgRXN0aW1hdGluZyB0aGUgQWNjdXJhY3kgb2YgYSBTdGF0aXN0aWMgb2YgSW50ZXJlc3QNClVuYSBkZSBsYXMgZ3JhbmRlcyB2ZW50YWphcyBkZWwgZW5mb3F1ZSBib290c3RyYXAgZXMgcXVlIHB1ZWRlIHNlciBhcGxpY2FkbyBlbiBjYXNpIHRvZGFzIGxhcyBzaXR1YWNpb25lcy4gTG9zIGNhbGN1bG9zIG1hdGVtYXRpY29zIGNvbXBsaWNhZG9zIG5vIHNvbiByZXF1ZXJpZG9zLiBSZWFsaXphciB1biBhbmFsaXNpcyBib290c3RyYXAgZW4gUiBjb25sbGV2YSBzb2xvIGRvcyBwYXNvcy4gDQoNCiAgLSBQcmltZXJvIHNlIGRlYmUgY3JlYXIgdW5hIGZ1bmNpb24gcXVlIGNvbXB1dGUgZWwgZXN0YWRpc3RpY28gZGUgaW50ZXJlcy4NCiAgLSBTZWd1bmRvIHVzYW1vcyBsYSBmdW5jaW9uIGBib290KClgIChkZSBsYSBsaWJyZXJpYSBgYm9vdGApIHBhcmEgcmVhbGl6YXIgZWwgYm9vdHN0cmFwLg0KDQpFbCBkYXRhc2V0IGBQb3J0YWZvbGlvYCAgc2UgdXRpbGl6YXJhIGVuIGVzdGUgZWplbXBsby4gUGFyYSBpbHVzdHJhciBlbCB1c28gZGUgYm9vdHN0cmFwIGVuIGVzdGUgZGF0YXNldCBwcmltZXJvIHNlIGRlYmUgY3JlYXIgdW5hIGZ1bmNpb24gYGFscGhhLmZuKClgIHF1ZSB0b21hICQoWCxZKSQgY29tbyBpbnB1dCBjb21vIHVuIHZlY3RvciBpbmRpY2FuZG8gY3VhbGVzIG9ic2VydmFjaW9uZXMgZGViZXJpYW4gc2VyIHVzYWRhcyBwYXJhIGVzdGltYXIgJFxhbHBoYSQuDQoNCmBgYHtyfQ0KYWxwaGEuZm4gPC0gZnVuY3Rpb24oZGF0YSxpbmRleCl7DQogIFggPC0gZGF0YSRYW2luZGV4XQ0KICBZIDwtIGRhdGEkWVtpbmRleF0NCiAgcmV0dXJuKCh2YXIoWSkgLSBjb3YoWCxZKSkgLyAodmFyKFgpICsgdmFyKFkpIC0gMipjb3YoWCxZKSkpDQp9DQpgYGANCg0KRXN0YSBmdW5jaW9uIGRldnVlbHZlIGNvbW8gb3V0cHV0IHVuIGVzdGltYWRvIHBhcmEgJFxhbHBoYSQgYmFzYWRvIGVuIGVsIGFyZ3VtZW50byBgaW5kZXhgLiBBaG9yYSBzZSB1c2FyYSBwYXJhIGVzdGltYXIgJFxhbHBoYSQgdXNhbmRvIDEwMCBvYnNlcnZhY2lvbmVzOg0KYGBge3J9DQphbHBoYS5mbihQb3J0Zm9saW8sIDE6MTAwKQ0KYGBgDQoNCkFob3JhIHNlIGhhcmEgdW4gYHNhbXBsZSgpYCBwYXJhIHNlbGVjY2lvbmFyIDEwMCBvYnNlcnZhY2lvbmVzIGRlbCByYW5nbyAxIGEgMTAwLiANCmBgYHtyfQ0Kc2V0LnNlZWQoMSkNCmFscGhhLmZuKFBvcnRmb2xpbywgc2FtcGxlKDEwMCwxMDAscmVwbGFjZSA9IFRSVUUpKQ0KYGBgDQoNClNlIHB1ZWRlIGltcGxlbWVudGFyIHVuIGFuYWxpc2lzIGJvb3RzdHJhcCByZWFsaXphbmRvIGVzdGUgY29tYW5kbyBtdWNoYXMgdmVjZXMuIEVzdG8gbG8gcHVlZGUgaGFjZXIgbGEgZnVuY2lvbiBgYm9vdCgpYCwgc2UgcHJvZHVjaXJhbiAkUj0xMDAwJCBlc3RpbWFjaW9uZXMgYm9vdHN0cmFwIHBhcmEgJFxhbHBoYSQuDQpgYGB7cn0NCmJvb3QoUG9ydGZvbGlvLCBhbHBoYS5mbiwgUiA9IDEwMDApDQpgYGANCg0KRWwgb3V0cHV0IG11ZXN0cmEgcXVlIHVzYW5kbyBlbCBkYXRhc2V0IG9yaWdpbmFsICRcaGF0e1xhbHBoYX0gPSAwLjU3NTgkIHkgcXVlIGVsIGJvb3RzdHJhcCBlc3RpbWFkbyBlcyAkU0UoXGhhdHtcYWxwaGF9PSAwLjA4ODYkLg0KDQojIyMgRXN0aW1hdGluZyB0aGUgQWNjdXJhY3kgb2YgYSBMaW5lYXIgUmVncmVzc2lvbiBNb2RlbA0KDQpFbCBlbmZvcXVlIGJvb3RzdHJhcCBwdWVkZSBzZXIgdXNhZG8gcGFyYSBldmFsdWFyIGxhIHZhcmlhYmxpZGFkIGRlIGxhcyBlc3RpbWFjaW9uZXMgZGUgbG9zIGNvZWZpY2llbnRlcyB5IGxhcyBwcmVkaWNjaW9uZXMgZGVsIG1ldG9kbyBkZSBhcHJlbmRpemFqZSBlc3RhZGlzdGljby4gDQpTZSB1dGlsaXphcmEgYm9vdHN0cmFwIHBhcmEgZXZhbHVhciBsYSB2YXJpYWJpbGlkYWQgZGUgbG9zIGVzdGltYWRvcyBwYXJhICRcYmV0YV8wJCB5ICRcYmV0YV8xJCwgZWwgaW50ZXJjZXB0byB5IGxhIHBlbmRpZW50ZSBwYXJhIGVsIG1vZGVsbyBkZSByZWdyZXNpb24gbGluZWFsIHF1ZSB1c2EgYGhvcnNlcG93ZXJgIHBhcmEgcHJlZGVjaXIgYG1wZ2AgZW4gZWwgZGF0YXNldCBgQXV0b2AuDQoNClByaW1lcm8gc2UgdmEgYSBjcmVhciB1bmEgZnVuY2lvbiBzaW1wcmUgYGJvb3QuZm4oKWAsIGxhIGN1YWwgdG9tYSBlbCBkYXRhc2V0IGBBdXRvYCB5IGxvcyBpbmRpY2VzIGRlIG9ic2VydmFjaW9uIHkgcmVncmVzYSBlbCBpbnRlcmNlcHRvIHkgbGEgcGVuZGllbnRlIGVzdGltYWRhcyBwYXJhIGVsIG1vZGVsbyBkZSByZWdyZXNpb24gbGluZWFsLiBEZXNwdWVzIGFwbGljYW1vcyBlc3RhIGZ1bmNpb24gYWwgZGF0YXNldCBjb21wbGV0byB5IGNvbXB1dGFyZW1vcyBsb3MgZXNtaW1hZG9zIGRlICRcYmV0YV8wJCB5IGRlICRcYmV0YV8xJC4NCmBgYHtyfQ0KYm9vdC5mbiA8LSBmdW5jdGlvbihkYXRhLCBpbmRleCkNCiAgcmV0dXJuKGNvZWYobG0obXBnIH4gaG9yc2Vwb3dlciwgZGF0YSA9IGRhdGEsIHN1YnNldCA9IGluZGV4KSkpDQpib290LmZuKEF1dG8sIDE6MzkyKQ0KYGBgDQoNCkxhIGZ1bmNpb24gYGJvb3QuZm4oKWAgcHVlZGUgc2VyIHVzYXJhIHBhcmEgY3JlYXIgZXN0aW1hZG9yZXMgYm9vdHN0cmFwIHBhcmEgbG9zIGludGVyY2VwdG9zIHkgbGEgcGVuZGllbnRlIGhhY2llbmRvIHVuIG11ZXN0cmVvIGRlIGxhcyBvYnNlcnZhY2lvbmVzOg0KYGBge3J9DQpzZXQuc2VlZCgxKQ0KYm9vdC5mbihBdXRvLCBzYW1wbGUoMzkyLCAzOTIsIHJlcGxhY2UgPSBUUlVFKSkNCmJvb3QuZm4oQXV0byxzYW1wbGUoMzkyLCAzOTIsIHJlcGxhY2UgPSBUUlVFKSkNCmBgYA0KDQpBaG9yYSB1c2FyZW1vcyBsYSBmdW5jaW9uIGBib290KClgIHBhcmEgb2J0ZW5lciBsb3MgZXJyb3JlcyBlc3RhbmRhciBkZSAkMTAwMCQgZXN0aW1hY2lvbmVzIGJvb3RzdHJhcCBwYXJhIGVsIGludGVyY2VwdG8geSBsYSBwZW5kaWVudGUuDQpgYGB7cn0NCmJvb3QoQXV0bywgYm9vdC5mbiwgMTAwMCkNCmBgYA0KDQpFc3RvIGluZGljYSBxdWUgZWwgZXN0aW1hZG9yIGJvb3RzdHJhcCBlcyAkU0UoXGhhdHtcYmV0YV8wfSkgPSAwLjg2JCB5ICRTRShcaGF0e1xiZXRhXzF9KSA9IDAuMDA3NCQuDQoNCkxhcyBmb3JtdWxhcyBlc3RhcmRhciBwdWVkZW4gc2VyIHVzYWRhcyBwYXJhIG9idGVuZXIgbG9zIGVycm9lcyBlc3RhbmRhciBkZSBsb3MgY29lZmljaWVudGVzIHJlZ3Jlc29yZXMgZW4gdW4gbW9kZWxvIGxpbmVhbDoNCmBgYHtyfQ0Kc3VtbWFyeShsbShtcGcgfiBob3JzZXBvd2VyLCBkYXRhID0gQXV0bykpJGNvZWYNCmBgYA0KDQpgYGB7cn0NCmJvb3QuZm4gPC0gZnVuY3Rpb24oZGF0YSwgaW5kZXgpDQogIGNvZWZmaWNpZW50cyhsbShtcGcgfiBob3JzZXBvd2VyICsgSShob3JzZXBvd2VyXjIpLCBkYXRhID0gZGF0YSwgc3Vic2V0ID0gaW5kZXgpKQ0Kc2V0LnNlZWQoMSkNCmJvb3QoQXV0bywgYm9vdC5mbiwgMTAwMCkNCnN1bW1hcnkobG0oIG1wZyB+IGhvcnNlcG93ZXIgKyBJKGhvcnNlcG93ZXJeMiksIGRhdGEgPSBBdXRvKSkkY29lZg0KYGBgDQoNCg0KDQoNCg0KDQoNCg==