Tutorials › Exam PA - Predictive Analytics › Worked Example: Fitting and Interpreting Ridge in R

Exam PA · Penalties and Ridge Regression · Tutorial 425 of 966

Worked Example: Fitting and Interpreting Ridge in R

Work through a reproducible glmnet ridge fit and learn how to carry its coefficients and predictions back to the units of the original predictors.

Intermediate 10 min read

What You'll Learn

  • Prepare a numeric model matrix and center and scale its predictors using training data.
  • Fit a Gaussian ridge model with glmnet and identify its penalty and intercept settings.
  • Extract and verify fitted coefficients against the ridge objective.
  • Convert standardized slopes and the intercept to original predictor units.
  • Make a prediction using consistent scaling and check it on both scales.
  • Avoid common errors involving intercepts, scaling, and model-matrix columns.

From a Data Frame to an Interpretable Ridge Fit

A ridge model is useful only if the fitted object can be connected correctly to the data used to fit it. In R, that means preparing the predictor columns, making the scaling convention explicit, fitting the penalized model, extracting its coefficients, and applying the same preparation to any new observations. The final interpretation should usually return to the units that a decision-maker recognizes, such as dollars, years, or counts.

This tutorial uses a small hypothetical data set so that the glmnet fit can be checked by hand. The example is deliberately transparent: its predictors are centered and uncorrelated, allowing us to verify the fitted coefficients from the ridge objective. Real data need not have this convenient structure, but the R workflow is the same. This material is aligned with the SOA Exam PA (Predictive Analytics) syllabus.

Definition: For a Gaussian ridge model with \(n\) observations, design matrix \(Z\), response vector \(y\), and an unpenalized intercept \(\beta_0\), the glmnet objective when \(\alpha=0\) is $$ \frac{1}{2n}\sum_{i=1}^n (y_i-\beta_0-z_i^\mathsf{T}\beta)^2+\frac{\lambda}{2}\sum_{j=1}^p\beta_j^2. $$ Here \(\lambda\geq 0\) controls the penalty, and \(Z\) contains the predictor columns supplied to the fit. The intercept is not penalized. The penalty applies to the slopes, so the meaning and scale of each slope depend on how the columns of \(Z\) were prepared.

Worked Example: Fit, Extract, Predict, and Interpret

Worked Example: A Reproducible Ridge Fit in R

A hypothetical analyst has four records with two risk measurements, \(x_1\) and \(x_2\), and a continuous response \(y\), measured in dollars. Fit a Gaussian ridge model with \(\lambda=1\), then predict for a record with \(x_1=55\) and \(x_2=110\). Use a training-set scaling convention in which each predictor is centered and divided by its root mean square deviation, calculated with divisor \(n\).

Record\(x_1\)\(x_2\)\(y\)
1408097
24012099
36080101
460120103

Step 1: Prepare the model matrix. With numeric predictors, the formula \(~x_1+x_2\) creates an intercept column and two predictor columns. glmnet fits its own intercept, so remove the model-matrix intercept column before supplying the predictors. The following R commands create that matrix and carry out the stated training-set scaling.

train <- data.frame(x1=c(40,40,60,60), x2=c(80,120,80,120), y=c(97,99,101,103))
Xraw <- model.matrix(~ x1 + x2, data=train)[, -1, drop=FALSE]
mu <- colMeans(Xraw)
Xcenter <- sweep(Xraw, 2, mu, "-")
s <- sqrt(colMeans(Xcenter^2))
Z <- sweep(Xcenter, 2, s, "/")

The predictor means are \(\mu_1=50\) and \(\mu_2=100\). The centered columns are \((-10,-10,10,10)\) and \((-20,20,-20,20)\). Their root mean square deviations are \(s_1=10\) and \(s_2=20\), respectively. Therefore the columns of \(Z\) are \((-1,-1,1,1)\) and \((-1,1,-1,1)\). This scaling uses divisor \(n=4\), as specified; it is not the sample-standard-deviation convention with divisor \(n-1\).

Step 2: Fit and extract the model. Set \(\alpha=0\) for ridge and disable glmnet's internal standardization because the columns have already been scaled. The intercept remains enabled.

fit <- glmnet::glmnet(Z, train$y, family="gaussian", alpha=0,
  lambda=1, standardize=FALSE, intercept=TRUE)
coef(fit)

The fitted coefficients are an intercept of \(100\), a coefficient of \(1\) for the first scaled predictor, and a coefficient of \(0.5\) for the second. Thus the fitted model on the scaled columns is \(\widehat y=100+1z_1+0.5z_2\).

Step 3: Check the fit. The response mean is \(100\), so the intercept is \(100\) because both predictor columns have mean zero. Also, \(Z^\mathsf{T}Z/n\) is the \(2\)-by-\(2\) identity matrix, and \(Z^\mathsf{T}(y-\bar y)/n=(2,1)^\mathsf{T}\). The ridge slope estimate is therefore \((I+\lambda I)^{-1}(2,1)^\mathsf{T}=(1,0.5)^\mathsf{T}\) when \(\lambda=1\). This agrees with the extracted glmnet coefficients.

Step 4: Predict using the fitted scale. For \(x_1=55\) and \(x_2=110\), the scaled values are \(z_1=(55-50)/10=0.5\) and \(z_2=(110-100)/20=0.5\). The prediction is \(100+(1)(0.5)+(0.5)(0.5)=100.75\) dollars.

In R, prepare the new observation with the training means and scales, not statistics calculated from the new observation:

newZ <- matrix(c((55-mu[1])/s[1], (110-mu[2])/s[2]), nrow=1)
predict(fit, newx=newZ, s=1)

This returns \(100.75\). The value of \(s\) in the predict call specifies the penalty parameter \(\lambda\); it is not either predictor's scaling value.

Step 5: Interpret on the original scale. The original-unit slope for \(x_j\) is the scaled slope divided by \(s_j\). Thus the slopes are \(1/10=0.1\) dollars per unit of \(x_1\) and \(0.5/20=0.025\) dollars per unit of \(x_2\). The original-scale intercept is \(100-(0.1)(50)-(0.025)(100)=92.5\). The equivalent equation is \(\widehat y=92.5+0.1x_1+0.025x_2\). Substitution gives \(92.5+0.1(55)+0.025(110)=92.5+5.5+2.75=100.75\), checking the prediction.

Method and mistake avoided. The method is to scale from training data, fit with internal scaling disabled, then reverse the scaling algebraically. A common mistake is to divide the slopes by their scale factors but leave the intercept at \(100\). That would give \(100+0.1x_1+0.025x_2\), which predicts \(108.25\) at the stated observation and is not the same model. The intercept must change when the predictors are returned to their original origins.

Why the Coefficient Conversion Works

Suppose a fitted model uses \(z_j=(x_j-\mu_j)/s_j\), where \(\mu_j\) and \(s_j>0\) are the training mean and scale for predictor \(j\). If the fitted scaled equation is \(\widehat y=b_0^*+\sum_j b_j^*z_j\), substitute the definition of \(z_j\):

$$ \widehat y=b_0^*+\sum_j b_j^*\frac{x_j-\mu_j}{s_j} =\left(b_0^*-\sum_j\frac{b_j^*\mu_j}{s_j}\right) +\sum_j\frac{b_j^*}{s_j}x_j. $$

Matching the coefficient of each \(x_j\) and the constant term gives \(b_j=b_j^*/s_j\) and \(b_0=b_0^*-\sum_j b_j\mu_j\). This is an algebraic change of coordinates, not a refit. It preserves every fitted value, provided the same means and scales used at fitting are used in the conversion.

The intercept conversion matters even when the fitted intercept has a simple interpretation on the centered scale. With centered predictors, it is the fitted response at the training means. On the original scale, it is the fitted response when all original predictors equal zero, which may be outside the observed range and may not have a useful practical interpretation.

Worked Example: Convert a Different Fitted Model

Worked Example: Return to Original Predictor Units

Suppose a ridge fit uses scaled predictors with training means \(\mu_1=40\), \(\mu_2=200\), and scales \(s_1=10\), \(s_2=50\). Its extracted scaled equation is \(\widehat y=12+3z_1-2z_2\). Convert it to original units and predict at \(x_1=50\), \(x_2=150\).

Method. Divide each scaled slope by its own scale: \(b_1=3/10=0.3\) and \(b_2=-2/50=-0.04\). Then convert the intercept: \(b_0=12-(0.3)(40)-(-0.04)(200)=12-12+8=8\). The original-scale equation is \(\widehat y=8+0.3x_1-0.04x_2\).

At the specified observation, \(z_1=(50-40)/10=1\) and \(z_2=(150-200)/50=-1\). The scaled equation gives \(12+3(1)-2(-1)=17\). The original equation gives \(8+0.3(50)-0.04(150)=8+15-6=17\). The two calculations agree.

Method and mistake avoided. Transform one slope at a time, then use the transformed slopes in the intercept formula. Merely subtracting the means from the old intercept, without multiplying by the converted slopes, would be incorrect.

Prediction and R Workflow Checks

The arithmetic in a hand calculation is only one part of a reliable R workflow. For a new record, the design matrix must have the same columns, in the same order, as the training matrix. When a predictor is a factor, model.matrix creates indicator columns according to the training factor levels and contrast settings. A new-data matrix must preserve those columns. A level absent from training data cannot be assigned a learned effect by the fitted model; do not silently treat such a level as an existing one.

For numeric predictors, the safest approach is to save the training means and scales and reuse them. Recomputing them on a test set or on an individual new record changes the predictor coordinates and therefore changes the meaning of the fitted coefficients. Similarly, if glmnet is allowed to standardize predictors internally, do not also assume that coefficients are on the scale you manually created. Make the choice explicit with the standardize argument.

Worked Example: Check a Prediction at Another Penalty

Consider centered, scaled predictor columns for which \(Z^\mathsf{T}Z/n=I\), and suppose \(Z^\mathsf{T}(y-\bar y)/n=(4,-2)^\mathsf{T}\) with \(\bar y=20\). Fit at \(\lambda=1\). The scales and means used to create the columns were \(s_1=2\), \(s_2=10\), \(\mu_1=10\), and \(\mu_2=100\). Predict at \(x_1=11\), \(x_2=90\).

Method. For orthogonal unit-RMS columns, the ridge slope equation is \((I+\lambda I)\widehat\beta=Z^\mathsf{T}(y-\bar y)/n\). Therefore the scaled slopes are \((4,-2)/(1+1)=(2,-1)\), and the intercept is \(20\). The new scaled values are \(z_1=(11-10)/2=0.5\) and \(z_2=(90-100)/10=-1\). The scaled prediction is \(20+2(0.5)+(-1)(-1)=22\).

Convert the slopes: \(b_1=2/2=1\) and \(b_2=-1/10=-0.1\). The original intercept is \(20-(1)(10)-(-0.1)(100)=20\). Thus the original equation is \(\widehat y=20+x_1-0.1x_2\), which gives \(20+11-9=22\). Both scales give the same prediction.

Method and mistake avoided. This check uses the specified penalty and then applies the saved training transformation to the new observation. A frequent error is to pass raw values \(11\) and \(90\) directly to a model fitted with scaled columns; that mixes coordinate systems and does not evaluate the fitted equation.

Practical Interpretation and Common Traps

A ridge coefficient describes a fitted association conditional on the other included predictors; it does not, by itself, establish a causal effect. On the scaled predictor scale, a one-unit increase means an increase of one training-set scale unit. On the original scale, the converted slope describes the fitted response change for a one-unit increase in that original predictor, holding the other model-matrix columns fixed. If predictors are correlated, that conditional interpretation can be sensitive to their joint information, even though ridge can stabilize the numerical fit.

Several implementation details deserve an explicit check before reporting results:

  • Intercept column: model.matrix normally includes a column of ones. Remove it when glmnet is also fitting an intercept; otherwise the model contains a redundant predictor column.
  • Penalty value: alpha selects the penalty family; alpha=0 specifies ridge. The lambda argument selects the penalty strength, and predict's s argument requests predictions at a penalty value.
  • Scaling choice: record the training means and scale factors, including whether the scale uses divisor \(n\) or \(n-1\). Coefficient conversion must match that choice.
  • Factor columns: retain the training model-matrix column definitions and their order when constructing new observations. The coefficient for an indicator is relative to its coding and reference level.
  • Units and interpretation: convert coefficients before describing effects in original units, and do not interpret the original-scale intercept as a typical response unless zero-valued predictors are meaningful.

The ridge objective provides a useful final check: with centered predictors, its intercept is the response mean, and its slopes balance fit to the response against the squared-coefficient penalty. For the special case of orthogonal columns with \(Z^\mathsf{T}Z/n=I\), the slope vector is \((1+\lambda)^{-1}Z^\mathsf{T}(y-\bar y)/n\). The calculation is not generally this simple for correlated columns, but extracted coefficients, predictions, and scale conversions can still be checked by applying the fitted equation consistently.

Key takeaway: Prepare the model matrix and scaling from the training data, fit ridge with a clearly specified lambda and standardization setting, and use the saved transformation for prediction. Convert both slopes and the intercept when reporting the model in original units.

Check Your Understanding

Work through each item, keeping the training transformation distinct from the prediction calculation.

  1. A scaled ridge equation is \(\widehat y=6+2z_1-z_2\), with training means \(30\) and \(40\), and scales \(5\) and \(10\). Give the original-scale intercept and slopes. (Answer: intercept \(-2\), slopes \(0.4\) and \(-0.1\).)
  2. For orthogonal scaled columns with \(Z^\mathsf{T}Z/n=I\), suppose the response-centered cross-product vector is \((3,-1)^\mathsf{T}\), the response mean is \(8\), and \(\lambda=2\). Give the fitted scaled slopes and intercept. (Answer: slopes \(1\) and \(-1/3\), intercept \(8\).)
  3. An original-scale fitted equation is \(\widehat y=15+0.2x_1+0.05x_2\). Calculate the prediction at \(x_1=20\), \(x_2=60\). (Answer: \(22\).)
  4. A model was trained with means \(10,20\), scales \(2,4\), and scaled equation \(\widehat y=5+3z_1+2z_2\). For a new record \(x_1=12,x_2=16\), calculate its scaled values and prediction. (Answer: \(z_1=1\), \(z_2=-1\), prediction \(6\).)
  5. For a centered, scaled fit, the extracted intercept is \(10\), the scaled slopes are \(2,-3\), and the training means and scales are \((5,20)\) and \((2,10)\), respectively. Convert the equation to original units. (Answer: intercept \(11\), slopes \(1,-0.3\), so \(\widehat y=11+x_1-0.3x_2\).)
Clinical Biostats