This post provides code to implement inductive selection/specification of log-linear models for contingency tables via the Lasso. I will exemplify the step-by-step implementation using data on intergenerational occupational mobility for 16 countries.

Steps:

  1. Install and load package glmnet, which aperforms Lasso or elastic-net regularization path for generalized linear models. Other packages for data manipulation will be also used.
    library("glmnet")
    library("tidyverse")
    library("modelr")
    library("reshape2")


  1. Input the contingency table an turn it into a data frame. I will use the dataset hg16 from the package logmult. This is a cross-classification of subject’s occupational status (destination) and his father’s occupational status (origin) across 16 countries.
    # Inpute data and create contingency table as data.frame()
    library("logmult")
    data(hg16)
    table <- ftable(hg16)
    mydata <- as.data.frame(table)
    names(mydata) <- c("origin","destination","country","Freq")

This is what the data looks like:


  1. Next, create the design matrix for the saturated sodel. In this example the saturated models contains 144 parameters.
    # Set reference categories
    mydata$origin      <- relevel(mydata$origin,ref = "Farm")
    mydata$destination <- relevel(mydata$destination,ref = "Farm")
    mydata$country     <- relevel(mydata$country,ref = "Spain")
    # Outcome variable is the frequencies in the contingency table
    y <- mydata$Freq + 0.5 # add small constant to avoid problems with empty cells. 
    # Creates design matrix with all three-way interactions
    formula <- as.formula( ~ .*.*.)
    x <- model.matrix(formula, mydata[,c("origin","destination","country")])


  1. Create a vector of penalty factors.
    # Penalty factors for weighted penalties. 
    # Here all variables have a penalty of 1, which is equivalent to unweighted penalties. 
    w <- rep(1,dim(x)[2])
    # Unsilence the lines below if you wish to implement the "adaptive lasso" based on ridge estimates
    #ridge.cv <- cv.glmnet(x,y,alpha=0,family='poisson', nfolds = 10)
    #best_ridge_coef <- as.numeric(coef(ridge.cv, s = ridge.cv$lambda.min))[-1] 
    #w <- 1/(abs(best_ridge_coef))^(1)


  1. Fit the saturated model using the glmnet package. In this case you should set the parameter \(\alpha=1\) so that you use a Lasso penalty (alternatively, \(\alpha=0\) corresponds to a ridge penalty and \(\alpha=0.5\) is elastic nets). Because the dependent variable is a vector of frequencies in a contingency table I use a Poisson family with a log-link function.
    # "x" is the design matrix for the saturated model 
    # "y" is the  dependent variable
    # "alpha=1" indicates the use of a Lasso penalty
    # "family=poisson" uses a log-link function by default
    # "penalty.factor = w specifies the weight for each variable's penalty
    
    # Fit model
    
    lasso <- glmnet(x,y,alpha=1,family='poisson', penalty.factor = w)
    # Plot the path of coefficients for values of lambda
    plot(lasso,xvar="lambda")


  1. Select a value of \(\lambda\) that minimizes a measure of cross-validation error. In the case of a Poisson model such measure is based on the Poisson Deviance. For this we use 10-fold cross-validation. Following the advice of the package developers I select the value of \(\lambda\) that yields the most regularized model such that the error is within one standard error of the minimum.
    # Cross-validation of Lambda
    lasso.cv <- cv.glmnet(x,y,alpha=1,family='poisson', nfolds = 10, penalty.factor = w)
    # Plot cross-validation error for each value of lambda
    plot(lasso.cv)


  1. Extract coefficients corresponding to the chosen value of lambda. The vector of coefficients from the Poisson model is most likely sparse.
# Extract optimum value of lambda 
opt.lam <- c(lasso.cv$lambda.1se)
lasso.coefs <- coef(lasso.cv, s = opt.lam)
print(lasso.coefs)
145 x 1 sparse Matrix of class "dgCMatrix"
                                                                          1
(Intercept)                                                      9.12085673
(Intercept)                                                      .         
originWhite Collar                                              -2.08988473
originBlue Collar                                               -1.91262871
destinationWhite Collar                                         -0.64397157
destinationBlue Collar                                          -0.13131034
countryAustralia                                                -3.08419718
countryBelgium                                                  -3.15994469
countryFrance                                                   -1.71517333
countryHungary                                                  -0.97210452
countryItaly                                                    -3.29253459
countryJapan                                                    -2.94361973
countryPhilippines                                              -1.17843790
countryUnited States                                            -2.03696358
countryWest Germany                                             -1.43042763
countryWest Malaysia                                            -1.47979006
countryYugoslavia                                               -3.78262426
countryDenmark                                                  -4.03057779
countryFinland                                                  -4.09436467
countryNorway                                                   -4.05957652
countrySweden                                                   -4.03704431
originWhite Collar:destinationWhite Collar                       2.49071711
originBlue Collar:destinationWhite Collar                        1.61331189
originWhite Collar:destinationBlue Collar                        0.73290796
originBlue Collar:destinationBlue Collar                         1.95783664
originWhite Collar:countryAustralia                              .         
originBlue Collar:countryAustralia                               0.07156901
originWhite Collar:countryBelgium                                .         
originBlue Collar:countryBelgium                                 .         
originWhite Collar:countryFrance                                 .         
originBlue Collar:countryFrance                                  .         
originWhite Collar:countryHungary                               -1.03005913
originBlue Collar:countryHungary                                 .         
originWhite Collar:countryItaly                                  .         
originBlue Collar:countryItaly                                   .         
originWhite Collar:countryJapan                                  .         
originBlue Collar:countryJapan                                   .         
originWhite Collar:countryPhilippines                           -0.54757703
originBlue Collar:countryPhilippines                            -0.69820992
originWhite Collar:countryUnited States                          .         
originBlue Collar:countryUnited States                           .         
originWhite Collar:countryWest Germany                           0.09309360
originBlue Collar:countryWest Germany                            .         
originWhite Collar:countryWest Malaysia                          .         
originBlue Collar:countryWest Malaysia                          -0.26259721
originWhite Collar:countryYugoslavia                             .         
originBlue Collar:countryYugoslavia                              .         
originWhite Collar:countryDenmark                                .         
originBlue Collar:countryDenmark                                 .         
originWhite Collar:countryFinland                                .         
originBlue Collar:countryFinland                                 .         
originWhite Collar:countryNorway                                 .         
originBlue Collar:countryNorway                                  .         
originWhite Collar:countrySweden                                 .         
originBlue Collar:countrySweden                                  .         
destinationWhite Collar:countryAustralia                         .         
destinationBlue Collar:countryAustralia                          .         
destinationWhite Collar:countryBelgium                           .         
destinationBlue Collar:countryBelgium                            .         
destinationWhite Collar:countryFrance                            .         
destinationBlue Collar:countryFrance                             .         
destinationWhite Collar:countryHungary                          -0.74394797
destinationBlue Collar:countryHungary                            .         
destinationWhite Collar:countryItaly                             .         
destinationBlue Collar:countryItaly                              .         
destinationWhite Collar:countryJapan                             .         
destinationBlue Collar:countryJapan                             -0.02763977
destinationWhite Collar:countryPhilippines                      -1.34939116
destinationBlue Collar:countryPhilippines                       -1.28669413
destinationWhite Collar:countryUnited States                     0.14914382
destinationBlue Collar:countryUnited States                      0.47475214
destinationWhite Collar:countryWest Germany                      .         
destinationBlue Collar:countryWest Germany                      -0.16409332
destinationWhite Collar:countryWest Malaysia                    -1.08710815
destinationBlue Collar:countryWest Malaysia                     -0.93576951
destinationWhite Collar:countryYugoslavia                        .         
destinationBlue Collar:countryYugoslavia                         .         
destinationWhite Collar:countryDenmark                           .         
destinationBlue Collar:countryDenmark                            .         
destinationWhite Collar:countryFinland                           .         
destinationBlue Collar:countryFinland                            .         
destinationWhite Collar:countryNorway                            .         
destinationBlue Collar:countryNorway                             .         
destinationWhite Collar:countrySweden                            .         
destinationBlue Collar:countrySweden                             .         
originWhite Collar:destinationWhite Collar:countryAustralia      .         
originBlue Collar:destinationWhite Collar:countryAustralia       .         
originWhite Collar:destinationBlue Collar:countryAustralia       .         
originBlue Collar:destinationBlue Collar:countryAustralia        0.16874354
originWhite Collar:destinationWhite Collar:countryBelgium        0.21522315
originBlue Collar:destinationWhite Collar:countryBelgium         0.17705495
originWhite Collar:destinationBlue Collar:countryBelgium         .         
originBlue Collar:destinationBlue Collar:countryBelgium          .         
originWhite Collar:destinationWhite Collar:countryFrance         0.42074875
originBlue Collar:destinationWhite Collar:countryFrance          0.24230301
originWhite Collar:destinationBlue Collar:countryFrance          0.91479074
originBlue Collar:destinationBlue Collar:countryFrance           0.39779339
originWhite Collar:destinationWhite Collar:countryHungary        .         
originBlue Collar:destinationWhite Collar:countryHungary         0.35107067
originWhite Collar:destinationBlue Collar:countryHungary         .         
originBlue Collar:destinationBlue Collar:countryHungary         -0.14851797
originWhite Collar:destinationWhite Collar:countryItaly          .         
originBlue Collar:destinationWhite Collar:countryItaly           .         
originWhite Collar:destinationBlue Collar:countryItaly           .         
originBlue Collar:destinationBlue Collar:countryItaly            .         
originWhite Collar:destinationWhite Collar:countryJapan          .         
originBlue Collar:destinationWhite Collar:countryJapan           .         
originWhite Collar:destinationBlue Collar:countryJapan           .         
originBlue Collar:destinationBlue Collar:countryJapan           -0.12605413
originWhite Collar:destinationWhite Collar:countryPhilippines    .         
originBlue Collar:destinationWhite Collar:countryPhilippines     .         
originWhite Collar:destinationBlue Collar:countryPhilippines     .         
originBlue Collar:destinationBlue Collar:countryPhilippines      .         
originWhite Collar:destinationWhite Collar:countryUnited States  0.34321015
originBlue Collar:destinationWhite Collar:countryUnited States   1.02222311
originWhite Collar:destinationBlue Collar:countryUnited States   0.18588765
originBlue Collar:destinationBlue Collar:countryUnited States    0.37988355
originWhite Collar:destinationWhite Collar:countryWest Germany   0.62405862
originBlue Collar:destinationWhite Collar:countryWest Germany    0.05653652
originWhite Collar:destinationBlue Collar:countryWest Germany    0.46213902
originBlue Collar:destinationBlue Collar:countryWest Germany     .         
originWhite Collar:destinationWhite Collar:countryWest Malaysia -0.04383566
originBlue Collar:destinationWhite Collar:countryWest Malaysia   .         
originWhite Collar:destinationBlue Collar:countryWest Malaysia   .         
originBlue Collar:destinationBlue Collar:countryWest Malaysia   -0.16241454
originWhite Collar:destinationWhite Collar:countryYugoslavia     .         
originBlue Collar:destinationWhite Collar:countryYugoslavia      .         
originWhite Collar:destinationBlue Collar:countryYugoslavia      .         
originBlue Collar:destinationBlue Collar:countryYugoslavia       .         
originWhite Collar:destinationWhite Collar:countryDenmark        .         
originBlue Collar:destinationWhite Collar:countryDenmark         .         
originWhite Collar:destinationBlue Collar:countryDenmark         .         
originBlue Collar:destinationBlue Collar:countryDenmark          .         
originWhite Collar:destinationWhite Collar:countryFinland        .         
originBlue Collar:destinationWhite Collar:countryFinland         .         
originWhite Collar:destinationBlue Collar:countryFinland         .         
originBlue Collar:destinationBlue Collar:countryFinland          .         
originWhite Collar:destinationWhite Collar:countryNorway         .         
originBlue Collar:destinationWhite Collar:countryNorway          .         
originWhite Collar:destinationBlue Collar:countryNorway          .         
originBlue Collar:destinationBlue Collar:countryNorway           .         
originWhite Collar:destinationWhite Collar:countrySweden         .         
originBlue Collar:destinationWhite Collar:countrySweden          .         
originWhite Collar:destinationBlue Collar:countrySweden          .         
originBlue Collar:destinationBlue Collar:countrySweden           .         


  1. You can plot the results to easly visualize the margin-free association between origin and destination.

LS0tCnRpdGxlOiAiTGFzc28gcmVndWxhcml6YXRpb24gZm9yIHNlbGVsZWN0aW9uL3NwZWNpZmljYXRpb24gb2YgbG9nLWxpbmVhciBtb2RlbHMiCm91dHB1dDoKICBodG1sX25vdGVib29rOiBkZWZhdWx0Ci0tLQoKCmBgYHtyIHNldHVwLCBpbmNsdWRlPUZBTFNFfQprbml0cjo6b3B0c19jaHVuayRzZXQoZWNobyA9IFRSVUUpCmBgYAoKClRoaXMgcG9zdCBwcm92aWRlcyBjb2RlIHRvIGltcGxlbWVudCBpbmR1Y3RpdmUgc2VsZWN0aW9uL3NwZWNpZmljYXRpb24gb2YgbG9nLWxpbmVhciBtb2RlbHMgZm9yIGNvbnRpbmdlbmN5IHRhYmxlcyB2aWEgdGhlIExhc3NvLiBJIHdpbGwgZXhlbXBsaWZ5IHRoZSBzdGVwLWJ5LXN0ZXAgaW1wbGVtZW50YXRpb24gdXNpbmcgZGF0YSBvbiBpbnRlcmdlbmVyYXRpb25hbCBvY2N1cGF0aW9uYWwgbW9iaWxpdHkgZm9yIDE2IGNvdW50cmllcy4KCiMjIyMgU3RlcHM6CgoxLiAgSW5zdGFsbCBhbmQgbG9hZCBwYWNrYWdlIGBgYGdsbW5ldGBgYCwgd2hpY2ggYXBlcmZvcm1zIExhc3NvIG9yIGVsYXN0aWMtbmV0IHJlZ3VsYXJpemF0aW9uIHBhdGggZm9yIGdlbmVyYWxpemVkIGxpbmVhciBtb2RlbHMuIE90aGVyIHBhY2thZ2VzIGZvciBkYXRhIG1hbmlwdWxhdGlvbiB3aWxsIGJlIGFsc28gdXNlZC4KCmBgYHtyLCBtZXNzYWdlPUZBTFNFfQogICAgbGlicmFyeSgiZ2xtbmV0IikKICAgIGxpYnJhcnkoInRpZHl2ZXJzZSIpCiAgICBsaWJyYXJ5KCJtb2RlbHIiKQogICAgbGlicmFyeSgicmVzaGFwZTIiKQpgYGAKCjxicj4KCjIuIElucHV0IHRoZSBjb250aW5nZW5jeSB0YWJsZSBhbiB0dXJuIGl0IGludG8gYSBkYXRhIGZyYW1lLiBJIHdpbGwgdXNlIHRoZSBkYXRhc2V0IGBgYGhnMTZgYGAgZnJvbSB0aGUgcGFja2FnZSBgYGBsb2dtdWx0YGBgLiBUaGlzIGlzIGEgY3Jvc3MtY2xhc3NpZmljYXRpb24gb2Ygc3ViamVjdCdzIG9jY3VwYXRpb25hbCBzdGF0dXMgKGRlc3RpbmF0aW9uKSBhbmQgaGlzIGZhdGhlcidzIG9jY3VwYXRpb25hbCBzdGF0dXMgKG9yaWdpbikgYWNyb3NzIDE2IGNvdW50cmllcy4gCgpgYGB7ciwgbWVzc2FnZT1GQUxTRX0KICAgICMgSW5wdXRlIGRhdGEgYW5kIGNyZWF0ZSBjb250aW5nZW5jeSB0YWJsZSBhcyBkYXRhLmZyYW1lKCkKICAgIGxpYnJhcnkoImxvZ211bHQiKQogICAgZGF0YShoZzE2KQogICAgdGFibGUgPC0gZnRhYmxlKGhnMTYpCiAgICBteWRhdGEgPC0gYXMuZGF0YS5mcmFtZSh0YWJsZSkKICAgIG5hbWVzKG15ZGF0YSkgPC0gYygib3JpZ2luIiwiZGVzdGluYXRpb24iLCJjb3VudHJ5IiwiRnJlcSIpCmBgYAoKVGhpcyBpcyB3aGF0IHRoZSBkYXRhIGxvb2tzIGxpa2U6CgpgYGB7ciwgZWNobz1GQUxTRX0KICAgIHByaW50KG15ZGF0YSAlPiUgYXNfdGliYmxlKCkpCmBgYAoKPGJyPgoKMy4gTmV4dCwgY3JlYXRlIHRoZSBkZXNpZ24gbWF0cml4IGZvciB0aGUgc2F0dXJhdGVkIHNvZGVsLiBJbiB0aGlzIGV4YW1wbGUgdGhlIHNhdHVyYXRlZCBtb2RlbHMgY29udGFpbnMgMTQ0IHBhcmFtZXRlcnMuIAoKYGBge3J9CiAgICAjIFNldCByZWZlcmVuY2UgY2F0ZWdvcmllcwoKICAgIG15ZGF0YSRvcmlnaW4gICAgICA8LSByZWxldmVsKG15ZGF0YSRvcmlnaW4scmVmID0gIkZhcm0iKQogICAgbXlkYXRhJGRlc3RpbmF0aW9uIDwtIHJlbGV2ZWwobXlkYXRhJGRlc3RpbmF0aW9uLHJlZiA9ICJGYXJtIikKICAgIG15ZGF0YSRjb3VudHJ5ICAgICA8LSByZWxldmVsKG15ZGF0YSRjb3VudHJ5LHJlZiA9ICJTcGFpbiIpCgogICAgIyBPdXRjb21lIHZhcmlhYmxlIGlzIHRoZSBmcmVxdWVuY2llcyBpbiB0aGUgY29udGluZ2VuY3kgdGFibGUKCiAgICB5IDwtIG15ZGF0YSRGcmVxICsgMC41ICMgYWRkIHNtYWxsIGNvbnN0YW50IHRvIGF2b2lkIHByb2JsZW1zIHdpdGggZW1wdHkgY2VsbHMuIAoKICAgICMgQ3JlYXRlcyBkZXNpZ24gbWF0cml4IHdpdGggYWxsIHRocmVlLXdheSBpbnRlcmFjdGlvbnMKCiAgICBmb3JtdWxhIDwtIGFzLmZvcm11bGEoIH4gLiouKi4pCiAgICB4IDwtIG1vZGVsLm1hdHJpeChmb3JtdWxhLCBteWRhdGFbLGMoIm9yaWdpbiIsImRlc3RpbmF0aW9uIiwiY291bnRyeSIpXSkKYGBgCgo8YnI+Cgo0LiBDcmVhdGUgYSB2ZWN0b3Igb2YgcGVuYWx0eSBmYWN0b3JzLgoKYGBge3J9CiAgICAjIFBlbmFsdHkgZmFjdG9ycyBmb3Igd2VpZ2h0ZWQgcGVuYWx0aWVzLiAKICAgICMgSGVyZSBhbGwgdmFyaWFibGVzIGhhdmUgYSBwZW5hbHR5IG9mIDEsIHdoaWNoIGlzIGVxdWl2YWxlbnQgdG8gdW53ZWlnaHRlZCBwZW5hbHRpZXMuIAoKICAgIHcgPC0gcmVwKDEsZGltKHgpWzJdKQoKICAgICMgVW5zaWxlbmNlIHRoZSBsaW5lcyBiZWxvdyBpZiB5b3Ugd2lzaCB0byBpbXBsZW1lbnQgdGhlICJhZGFwdGl2ZSBsYXNzbyIgYmFzZWQgb24gcmlkZ2UgZXN0aW1hdGVzCgogICAgI3JpZGdlLmN2IDwtIGN2LmdsbW5ldCh4LHksYWxwaGE9MCxmYW1pbHk9J3BvaXNzb24nLCBuZm9sZHMgPSAxMCkKICAgICNiZXN0X3JpZGdlX2NvZWYgPC0gYXMubnVtZXJpYyhjb2VmKHJpZGdlLmN2LCBzID0gcmlkZ2UuY3YkbGFtYmRhLm1pbikpWy0xXSAKICAgICN3IDwtIDEvKGFicyhiZXN0X3JpZGdlX2NvZWYpKV4oMSkKCmBgYAoKPGJyPgoKNS4gRml0IHRoZSBzYXR1cmF0ZWQgbW9kZWwgdXNpbmcgdGhlIGBgYGdsbW5ldGBgYCBwYWNrYWdlLiBJbiB0aGlzIGNhc2UgeW91IHNob3VsZCBzZXQgdGhlIHBhcmFtZXRlciAkXGFscGhhPTEkIHNvIHRoYXQgeW91IHVzZSBhIExhc3NvIHBlbmFsdHkgKGFsdGVybmF0aXZlbHksICRcYWxwaGE9MCQgY29ycmVzcG9uZHMgdG8gYSByaWRnZSBwZW5hbHR5IGFuZCAkXGFscGhhPTAuNSQgaXMgZWxhc3RpYyBuZXRzKS4gQmVjYXVzZSB0aGUgZGVwZW5kZW50IHZhcmlhYmxlIGlzIGEgdmVjdG9yIG9mIGZyZXF1ZW5jaWVzIGluIGEgY29udGluZ2VuY3kgdGFibGUgSSB1c2UgYSBQb2lzc29uIGZhbWlseSB3aXRoIGEgbG9nLWxpbmsgZnVuY3Rpb24uIAogICAgCmBgYHtyfQogICAgIyAieCIgaXMgdGhlIGRlc2lnbiBtYXRyaXggZm9yIHRoZSBzYXR1cmF0ZWQgbW9kZWwgCiAgICAjICJ5IiBpcyB0aGUgIGRlcGVuZGVudCB2YXJpYWJsZQogICAgIyAiYWxwaGE9MSIgaW5kaWNhdGVzIHRoZSB1c2Ugb2YgYSBMYXNzbyBwZW5hbHR5CiAgICAjICJmYW1pbHk9cG9pc3NvbiIgdXNlcyBhIGxvZy1saW5rIGZ1bmN0aW9uIGJ5IGRlZmF1bHQKICAgICMgInBlbmFsdHkuZmFjdG9yID0gdyBzcGVjaWZpZXMgdGhlIHdlaWdodCBmb3IgZWFjaCB2YXJpYWJsZSdzIHBlbmFsdHkKICAgIAogICAgIyBGaXQgbW9kZWwKICAgIAogICAgbGFzc28gPC0gZ2xtbmV0KHgseSxhbHBoYT0xLGZhbWlseT0ncG9pc3NvbicsIHBlbmFsdHkuZmFjdG9yID0gdykKCiAgICAjIFBsb3QgdGhlIHBhdGggb2YgY29lZmZpY2llbnRzIGZvciB2YWx1ZXMgb2YgbGFtYmRhCgogICAgcGxvdChsYXNzbyx4dmFyPSJsYW1iZGEiKQpgYGAKCjxicj4KCjYuIFNlbGVjdCBhIHZhbHVlIG9mICRcbGFtYmRhJCB0aGF0IG1pbmltaXplcyBhIG1lYXN1cmUgb2YgY3Jvc3MtdmFsaWRhdGlvbiBlcnJvci4gSW4gdGhlIGNhc2Ugb2YgYSBQb2lzc29uIG1vZGVsIHN1Y2ggbWVhc3VyZSBpcyBiYXNlZCBvbiB0aGUgUG9pc3NvbiBEZXZpYW5jZS4gRm9yIHRoaXMgd2UgdXNlIDEwLWZvbGQgY3Jvc3MtdmFsaWRhdGlvbi4gRm9sbG93aW5nIHRoZSBhZHZpY2Ugb2YgdGhlIHBhY2thZ2UgZGV2ZWxvcGVycyBJIHNlbGVjdCB0aGUgdmFsdWUgb2YgJFxsYW1iZGEkIHRoYXQgeWllbGRzIHRoZSBtb3N0IHJlZ3VsYXJpemVkIG1vZGVsIHN1Y2ggdGhhdCB0aGUgZXJyb3IgaXMgd2l0aGluIG9uZSBzdGFuZGFyZCBlcnJvciBvZiB0aGUgbWluaW11bS4KCgpgYGB7cn0KCiAgICAjIENyb3NzLXZhbGlkYXRpb24gb2YgTGFtYmRhCgogICAgbGFzc28uY3YgPC0gY3YuZ2xtbmV0KHgseSxhbHBoYT0xLGZhbWlseT0ncG9pc3NvbicsIG5mb2xkcyA9IDEwLCBwZW5hbHR5LmZhY3RvciA9IHcpCgogICAgIyBQbG90IGNyb3NzLXZhbGlkYXRpb24gZXJyb3IgZm9yIGVhY2ggdmFsdWUgb2YgbGFtYmRhCgogICAgcGxvdChsYXNzby5jdikKYGBgCgo8YnI+Cgo3LiBFeHRyYWN0IGNvZWZmaWNpZW50cyBjb3JyZXNwb25kaW5nIHRvIHRoZSBjaG9zZW4gdmFsdWUgb2YgbGFtYmRhLiBUaGUgdmVjdG9yIG9mIGNvZWZmaWNpZW50cyBmcm9tIHRoZSBQb2lzc29uIG1vZGVsIGlzIG1vc3QgbGlrZWx5IHNwYXJzZS4gCiAKYGBge3J9CiMgRXh0cmFjdCBvcHRpbXVtIHZhbHVlIG9mIGxhbWJkYSAKCm9wdC5sYW0gPC0gYyhsYXNzby5jdiRsYW1iZGEuMXNlKQpsYXNzby5jb2VmcyA8LSBjb2VmKGxhc3NvLmN2LCBzID0gb3B0LmxhbSkKcHJpbnQobGFzc28uY29lZnMpCmBgYAoKPGJyPgoKOC4gWW91IGNhbiBwbG90IHRoZSByZXN1bHRzIHRvIGVhc2x5IHZpc3VhbGl6ZSB0aGUgbWFyZ2luLWZyZWUgYXNzb2NpYXRpb24gYmV0d2VlbiBvcmlnaW4gYW5kIGRlc3RpbmF0aW9uLiAKCgpgYGB7ciwgZWNobz1GQUxTRX0KCiMgQ3JlYXRlIHByZWRpY3Rpb25zIGZvciBtYXJnaW5zIGFuZCBMT1JhdGl2ZSBtYXRpbmcuIFVzZWQgaW4gcGxvdCAKCmR1bW15Lm1vZGVsIDwtIGxtKEZyZXEgfiBvcmlnaW4gKyBkZXN0aW5hdGlvbiArIGNvdW50cnksIGRhdGE9bXlkYXRhKQpuZXdfeCA8LSBteWRhdGEgJT4lIGRhdGFfZ3JpZChvcmlnaW4sZGVzdGluYXRpb24sY291bnRyeSwubW9kZWw9ZHVtbXkubW9kZWwpICU+JSBtb2RlbC5tYXRyaXgoZm9ybXVsYSwgLikKCgojIGZ1bGwgcHJlZGljdGlvbgpwcmVkaWN0aW9ucyA8LSBjYmluZChteWRhdGElPiUgZGF0YV9ncmlkKG9yaWdpbixkZXN0aW5hdGlvbixjb3VudHJ5LC5tb2RlbD1kdW1teS5tb2RlbCksIHByZWRpY3QobGFzc28uY3YsIG5ld194LCBzPW9wdC5sYW0pKSAlPiUKICBhc190aWJibGUoKSAlPiUgcmVuYW1lKHByZWQgPSBgMWApIAoKCiMgSW50ZXJjZXB0CmludGVyY2VwdCA8LSBwcmVkaWN0aW9ucyAlPiUgZmlsdGVyKG9yaWdpbj09IkZhcm0iLCBkZXN0aW5hdGlvbj09IkZhcm0iLCBjb3VudHJ5PT0iU3BhaW4iKSAlPiUgc3VtbWFyaXNlKHByZWQpICU+JSBhcy5udW1lcmljKCkKCgojIG1hcmdpbnMKcHJlZGljdGlvbnMgPC0gcHJlZGljdGlvbnMgJT4lIG11dGF0ZShwcmVkID0gcHJlZCAtIGludGVyY2VwdCkgIyByZW1vdmUgaW50ZXJjZXB0CgpwcmVkaWN0aW9uc19jb3VudHJ5ICAgICA8LSBwcmVkaWN0aW9ucyAlPiUgZmlsdGVyKG9yaWdpbj09IkZhcm0iLCBkZXN0aW5hdGlvbj09IkZhcm0iKSAgJT4lIHJlbmFtZShtYXJnaW5fY291bnRyeT1wcmVkKSAlPiUgc2VsZWN0KGNvdW50cnksbWFyZ2luX2NvdW50cnkpCnByZWRpY3Rpb25zX29yaWdpbiAgICAgIDwtIHByZWRpY3Rpb25zICU+JSBmaWx0ZXIoY291bnRyeT09IlNwYWluIiwgZGVzdGluYXRpb249PSJGYXJtIikgJT4lIHJlbmFtZShtYXJnaW5fb3JpZ2luPXByZWQpICU+JSBzZWxlY3Qob3JpZ2luLG1hcmdpbl9vcmlnaW4pCnByZWRpY3Rpb25zX2Rlc3RpbmF0aW9uIDwtIHByZWRpY3Rpb25zICU+JSBmaWx0ZXIoY291bnRyeT09IlNwYWluIiwgb3JpZ2luPT0iRmFybSIpICU+JSByZW5hbWUobWFyZ2luX2Rlc3RpbmF0aW9uPXByZWQpICU+JSBzZWxlY3QoZGVzdGluYXRpb24sbWFyZ2luX2Rlc3RpbmF0aW9uKQoKCiMgbWF0Y2gKcHJlZGljdGlvbnMgPC0gcHJlZGljdGlvbnMgJT4lIGxlZnRfam9pbihwcmVkaWN0aW9uc19jb3VudHJ5LCBieT0iY291bnRyeSIpCnByZWRpY3Rpb25zIDwtIHByZWRpY3Rpb25zICU+JSBsZWZ0X2pvaW4ocHJlZGljdGlvbnNfb3JpZ2luLCBieT0ib3JpZ2luIikKcHJlZGljdGlvbnMgPC0gcHJlZGljdGlvbnMgJT4lIGxlZnRfam9pbihwcmVkaWN0aW9uc19kZXN0aW5hdGlvbiwgYnk9ImRlc3RpbmF0aW9uIikKCgojIG1hcmdpbnMgYnkgY291bnRyeQpwcmVkaWN0aW9uc19jb3VudHJ5X29yaWdpbiA8LSBwcmVkaWN0aW9ucyAlPiUgZmlsdGVyKG9yaWdpbiE9IkZhcm0iLGNvdW50cnkhPSJTcGFpbiIsZGVzdGluYXRpb249PSJGYXJtIikgJT4lCiAgcmVuYW1lKG1hcmdpbl9jb3VudHJ5X29yaWdpbj1wcmVkKSAlPiUgbXV0YXRlKG1hcmdpbl9jb3VudHJ5X29yaWdpbiA9IG1hcmdpbl9jb3VudHJ5X29yaWdpbiAtIChtYXJnaW5fY291bnRyeSArIG1hcmdpbl9vcmlnaW4gKSkgJT4lCiAgc2VsZWN0KGNvdW50cnksb3JpZ2luLG1hcmdpbl9jb3VudHJ5X29yaWdpbikKCnByZWRpY3Rpb25zX2NvdW50cnlfZGVzdGluYXRpb24gPC0gcHJlZGljdGlvbnMgJT4lIGZpbHRlcihvcmlnaW49PSJGYXJtIixjb3VudHJ5IT0iU3BhaW4iLGRlc3RpbmF0aW9uIT0iRmFybSIpICU+JQogIHJlbmFtZShtYXJnaW5fY291bnRyeV9kZXN0aW5hdGlvbj1wcmVkKSAlPiUgbXV0YXRlKG1hcmdpbl9jb3VudHJ5X2Rlc3RpbmF0aW9uID0gbWFyZ2luX2NvdW50cnlfZGVzdGluYXRpb24gLSAobWFyZ2luX2NvdW50cnkgKyBtYXJnaW5fZGVzdGluYXRpb24gKSkgJT4lIAogIHNlbGVjdChjb3VudHJ5LGRlc3RpbmF0aW9uLG1hcmdpbl9jb3VudHJ5X2Rlc3RpbmF0aW9uKQoKcHJlZGljdGlvbnMgPC0gcHJlZGljdGlvbnMgJT4lIGxlZnRfam9pbihwcmVkaWN0aW9uc19jb3VudHJ5X29yaWdpbiwgYnk9YygiY291bnRyeSIsIm9yaWdpbiIpKSAlPiUgIHJlcGxhY2VfbmEobGlzdChtYXJnaW5fY291bnRyeV9vcmlnaW4gPSAwKSkKcHJlZGljdGlvbnMgPC0gcHJlZGljdGlvbnMgJT4lIGxlZnRfam9pbihwcmVkaWN0aW9uc19jb3VudHJ5X2Rlc3RpbmF0aW9uLCBieT1jKCJjb3VudHJ5IiwiZGVzdGluYXRpb24iKSkgJT4lICByZXBsYWNlX25hKGxpc3QobWFyZ2luX2NvdW50cnlfZGVzdGluYXRpb24gPSAwKSkKCgojIGNvbXB1dGVzIG1hcmdpbi1mcmVlIGxvZy1vZGQgcmF0aW9ucyAoTE9ScykKcHJlZGljdGlvbnMgPC0gcHJlZGljdGlvbnMgJT4lIAogIG11dGF0ZShMT1IgPSBwcmVkIC0gKG1hcmdpbl9jb3VudHJ5ICsgbWFyZ2luX29yaWdpbiArIG1hcmdpbl9kZXN0aW5hdGlvbiArIG1hcmdpbl9jb3VudHJ5X29yaWdpbiArIG1hcmdpbl9jb3VudHJ5X2Rlc3RpbmF0aW9uKSApIApsZXZlbHMub3JpZ2luICAgICAgPC0gYygiRmFybSIsIkJsdWUgQ29sbGFyIiwiV2hpdGUgQ29sbGFyIikKbGV2ZWxzLmRlc3RpbmF0aW9uIDwtIGMoIkZhcm0iLCJCbHVlIENvbGxhciIsIldoaXRlIENvbGxhciIpCmxldmVscy5jb3VudHJ5ICAgICA8LSBsZXZlbHMobXlkYXRhJGNvdW50cnkpCgoKcGxvdCA8LSBwcmVkaWN0aW9ucyAlPiUgCiAgZ2dwbG90KGFlcyh5PWZhY3RvcihvcmlnaW4sIGxldmVscyA9IHJldihsZXZlbHMub3JpZ2luKSksCiAgICAgICAgICAgICB4PWZhY3RvcihkZXN0aW5hdGlvbiwgbGV2ZWxzID0gbGV2ZWxzLmRlc3RpbmF0aW9uKSkpICsgZmFjZXRfd3JhcCggfiBjb3VudHJ5KSArIGdlb21fcmFzdGVyKGFlcyhmaWxsPSBMT1IpKSArCiAgc2NhbGVfZmlsbF9ncmFkaWVudG4obGltaXRzPWMoLTMuNSwzLjUpLCBjb2xvdXJzPWMoInJlZCIsIndoaXRlIiwiYmx1ZSIpKSArCiAgbGFicyh5PSJGYXRoZXIncyBvY2N1cGF0aW9uIiwgeD0gIkNoaWxkcmVuJ3Mgb2NjdXBhdGlvbiIsIGNvbG91cj0iIikgKwogIHRoZW1lX2J3KCkgKyB0aGVtZShheGlzLnRleHQueCA9IGVsZW1lbnRfdGV4dChzaXplPTksIGFuZ2xlPTQ1LCB2anVzdD0tMSwgaGp1c3Q9MCksCiAgICAgICAgICAgICAgICAgICAgIGF4aXMudGV4dC55ID0gZWxlbWVudF90ZXh0KHNpemU9OSwgYW5nbGU9MCksCiAgICAgICAgICAgICAgICAgICAgIHBsb3QudGl0bGU9IGVsZW1lbnRfdGV4dChzaXplPTExKSkgKwogIHNjYWxlX3hfZGlzY3JldGUocG9zaXRpb249InRvcCIpCgoKcHJpbnQocGxvdCkKCgpgYGAKCgo=