Penalized Regression on RSS (Summary Statistics) Objective
Source:R/regularizedRegressionWrappers.R
penalizedRss.RdGeneralizes lassosumRss() to support LASSO, MCP, SCAD, L0, L0L1,
and L0L2 penalties. Uses coordinate descent on the objective
$$\beta^T R \beta - 2 \beta^T z + \mathrm{penalty}(\beta)$$
where \(R\) is a (possibly pre-shrunk) LD matrix and \(z = \hat\beta /
\sqrt{n}\).
Arguments
- bhat
Numeric vector of marginal effect estimates (length p).
- R
The LD correlation matrix (a single matrix over the analysed window), as in
lassosumRss().- n
GWAS sample size (positive scalar).
- penalty
Penalty type:
"lasso","MCP","SCAD","L0","L0L1", or"L0L2".- lambda
Numeric vector of regularization parameter values along which to trace a solution path (warm-started, largest-first). For LASSO/MCP/SCAD this is the primary penalty strength; for L0 variants it controls the L1 component.
- gamma
Concavity parameter for MCP (default 3) or SCAD (default 3.7). Ignored for LASSO and L0 variants.
- alpha
Elastic-net mixing for MCP/SCAD: \(l_1 = \lambda \alpha\), \(l_2 = \lambda (1-\alpha)\). Default 1 (pure L1, no ridge).
- lambda0
L0 penalty weight (number of non-zeros). Required for L0 variants; ignored otherwise. Default 0.
- lambda2
L2 penalty weight for L0L2 variant. Default 0.
- thr
Convergence threshold. Default 1e-4.
- maxiter
Maximum coordinate descent iterations per lambda. Default 10000.
- maxSwaps
Maximum swap rounds for L0 variants. Default 100. Set to 0 to disable swaps.
Value
A list with components:
- beta
p x length(lambda) matrix of coefficient estimates.
- lambda
The lambda values used.
- conv
Convergence indicators (1 = converged).
- loss
Quadratic loss at each lambda.
- fbeta
Full penalized objective at each lambda.
- nparams
Number of non-zero coefficients at each lambda.
- betaEst
Coefficient vector at the lambda minimizing fbeta.
Examples
set.seed(42)
p <- 10; n <- 100
bhat <- rnorm(p, sd = 0.1)
R <- diag(p)
# MCP
penalizedRss(bhat, R, n, penalty = "MCP")
#> $beta
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 0.0137095845 0.0137095845 0.0137095845 0.0137095845 0.0137095845
#> [2,] -0.0056469817 -0.0056469817 -0.0056469817 -0.0056469817 -0.0056469817
#> [3,] 0.0036312841 0.0036312841 0.0036312841 0.0036312841 0.0036312841
#> [4,] 0.0063286260 0.0063286260 0.0063286260 0.0063286260 0.0063286260
#> [5,] 0.0040426832 0.0040426832 0.0040426832 0.0040426832 0.0040426832
#> [6,] -0.0010612452 -0.0010612452 -0.0010612452 -0.0010612452 -0.0009496679
#> [7,] 0.0151152200 0.0151152200 0.0151152200 0.0151152200 0.0151152200
#> [8,] -0.0009465904 -0.0009465904 -0.0009465904 -0.0009465904 -0.0007776857
#> [9,] 0.0201842371 0.0201842371 0.0201842371 0.0201842371 0.0201842371
#> [10,] -0.0006271410 -0.0006271410 -0.0006271410 -0.0004942588 -0.0002985116
#> [,6] [,7] [,8] [,9] [,10]
#> [1,] 1.370958e-02 1.370958e-02 0.013709584 0.013709584 0.013709584
#> [2,] -5.646982e-03 -5.646982e-03 -0.005646982 -0.005646982 -0.004515496
#> [3,] 3.631284e-03 3.631284e-03 0.003535514 0.002697455 0.001491950
#> [4,] 6.328626e-03 6.328626e-03 0.006328626 0.006328626 0.005537963
#> [5,] 4.042683e-03 4.042683e-03 0.004042683 0.003314554 0.002109048
#> [6,] -6.680954e-04 -2.630676e-04 0.000000000 0.000000000 0.000000000
#> [7,] 1.511522e-02 1.511522e-02 0.015115220 0.015115220 0.015115220
#> [8,] -4.961133e-04 -9.108539e-05 0.000000000 0.000000000 0.000000000
#> [9,] 2.018424e-02 2.018424e-02 0.020184237 0.020184237 0.020184237
#> [10,] -1.693917e-05 0.000000e+00 0.000000000 0.000000000 0.000000000
#> [,11] [,12] [,13] [,14] [,15] [,16]
#> [1,] 0.0137095845 0.0123809845 0.008792977 0.003631808 0.000000000 0
#> [2,] -0.0027814373 -0.0002870804 0.000000000 0.000000000 0.000000000 0
#> [3,] 0.0000000000 0.0000000000 0.000000000 0.000000000 0.000000000 0
#> [4,] 0.0038039038 0.0013095469 0.000000000 0.000000000 0.000000000 0
#> [5,] 0.0003749896 0.0000000000 0.000000000 0.000000000 0.000000000 0
#> [6,] 0.0000000000 0.0000000000 0.000000000 0.000000000 0.000000000 0
#> [7,] 0.0151152200 0.0144894378 0.010901430 0.005740262 0.000000000 0
#> [8,] 0.0000000000 0.0000000000 0.000000000 0.000000000 0.000000000 0
#> [9,] 0.0201842371 0.0201842371 0.018504956 0.013343787 0.005919705 0
#> [10,] 0.0000000000 0.0000000000 0.000000000 0.000000000 0.000000000 0
#> [,17] [,18] [,19] [,20]
#> [1,] 0 0 0 0
#> [2,] 0 0 0 0
#> [3,] 0 0 0 0
#> [4,] 0 0 0 0
#> [5,] 0 0 0 0
#> [6,] 0 0 0 0
#> [7,] 0 0 0 0
#> [8,] 0 0 0 0
#> [9,] 0 0 0 0
#> [10,] 0 0 0 0
#>
#> $lambda
#> [1] 0.0001000000 0.0001438450 0.0002069138 0.0002976351 0.0004281332
#> [6] 0.0006158482 0.0008858668 0.0012742750 0.0018329807 0.0026366509
#> [11] 0.0037926902 0.0054555948 0.0078475997 0.0112883789 0.0162377674
#> [16] 0.0233572147 0.0335981829 0.0483293024 0.0695192796 0.1000000000
#>
#> $conv
#> [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
#>
#> $loss
#> [1] -0.0009277110 -0.0009277110 -0.0009277110 -0.0009276934 -0.0009275620
#> [6] -0.0009269812 -0.0009259487 -0.0009252863 -0.0009238932 -0.0009150743
#> [11] -0.0008840717 -0.0008396894 -0.0007790770 -0.0005875828 -0.0002039265
#> [16] 0.0000000000 0.0000000000 0.0000000000 0.0000000000 0.0000000000
#>
#> $fbeta
#> [1] -9.274110e-04 -9.270903e-04 -9.264266e-04 -9.250887e-04 -9.225100e-04
#> [6] -9.177926e-04 -9.088671e-04 -8.912145e-04 -8.575437e-04 -7.997471e-04
#> [11] -7.092633e-04 -5.614662e-04 -3.590607e-04 -1.494649e-04 -2.336194e-05
#> [16] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
#>
#> $nparams
#> [1] 10 10 10 10 10 10 9 7 7 7 6 5 3 3 1 0 0 0 0 0
#>
#> $betaEst
#> [1] 0.0137095845 -0.0056469817 0.0036312841 0.0063286260 0.0040426832
#> [6] -0.0010612452 0.0151152200 -0.0009465904 0.0201842371 -0.0006271410
#>
# SCAD
penalizedRss(bhat, R, n, penalty = "SCAD")
#> $beta
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 0.0137095845 0.0137095845 0.0137095845 0.0137095845 0.0137095845
#> [2,] -0.0056469817 -0.0056469817 -0.0056469817 -0.0056469817 -0.0056469817
#> [3,] 0.0036312841 0.0036312841 0.0036312841 0.0036312841 0.0036312841
#> [4,] 0.0063286260 0.0063286260 0.0063286260 0.0063286260 0.0063286260
#> [5,] 0.0040426832 0.0040426832 0.0040426832 0.0040426832 0.0040426832
#> [6,] -0.0010612452 -0.0010612452 -0.0010612452 -0.0010377129 -0.0007536876
#> [7,] 0.0151152200 0.0151152200 0.0151152200 0.0151152200 0.0151152200
#> [8,] -0.0009465904 -0.0009465904 -0.0009465904 -0.0008556141 -0.0005715889
#> [9,] 0.0201842371 0.0201842371 0.0201842371 0.0201842371 0.0201842371
#> [10,] -0.0006271410 -0.0006271410 -0.0005457056 -0.0003482533 -0.0001990078
#> [,6] [,7] [,8] [,9] [,10]
#> [1,] 1.370958e-02 1.370958e-02 0.013709584 0.013709584 0.0137095845
#> [2,] -5.646982e-03 -5.646982e-03 -0.005646982 -0.004979307 -0.0032301425
#> [3,] 3.631284e-03 3.631284e-03 0.002993912 0.001798303 0.0009946332
#> [4,] 6.328626e-03 6.328626e-03 0.006328626 0.006061919 0.0043127541
#> [5,] 4.042683e-03 4.042683e-03 0.003647310 0.002431304 0.0014060323
#> [6,] -4.453969e-04 -1.753784e-04 0.000000000 0.000000000 0.0000000000
#> [7,] 1.511522e-02 1.511522e-02 0.015115220 0.015115220 0.0151152200
#> [8,] -3.307422e-04 -6.072359e-05 0.000000000 0.000000000 0.0000000000
#> [9,] 2.018424e-02 2.018424e-02 0.020184237 0.020184237 0.0201842371
#> [10,] -1.129278e-05 0.000000e+00 0.000000000 0.000000000 0.0000000000
#> [,11] [,12] [,13] [,14] [,15] [,16] [,17]
#> [1,] 0.013519367 0.0099001043 0.005861985 0.002421206 0.00000000 0 0
#> [2,] -0.001854292 -0.0001913869 0.000000000 0.000000000 0.00000000 0 0
#> [3,] 0.000000000 0.0000000000 0.000000000 0.000000000 0.00000000 0 0
#> [4,] 0.002535936 0.0008730313 0.000000000 0.000000000 0.00000000 0 0
#> [5,] 0.000249993 0.0000000000 0.000000000 0.000000000 0.00000000 0 0
#> [6,] 0.000000000 0.0000000000 0.000000000 0.000000000 0.00000000 0 0
#> [7,] 0.015115220 0.0121325843 0.007267620 0.003826841 0.00000000 0 0
#> [8,] 0.000000000 0.0000000000 0.000000000 0.000000000 0.00000000 0 0
#> [9,] 0.020184237 0.0201833762 0.014977248 0.008895858 0.00394647 0 0
#> [10,] 0.000000000 0.0000000000 0.000000000 0.000000000 0.00000000 0 0
#> [,18] [,19] [,20]
#> [1,] 0 0 0
#> [2,] 0 0 0
#> [3,] 0 0 0
#> [4,] 0 0 0
#> [5,] 0 0 0
#> [6,] 0 0 0
#> [7,] 0 0 0
#> [8,] 0 0 0
#> [9,] 0 0 0
#> [10,] 0 0 0
#>
#> $lambda
#> [1] 0.0001000000 0.0001438450 0.0002069138 0.0002976351 0.0004281332
#> [6] 0.0006158482 0.0008858668 0.0012742750 0.0018329807 0.0026366509
#> [11] 0.0037926902 0.0054555948 0.0078475997 0.0112883789 0.0162377674
#> [16] 0.0233572147 0.0335981829 0.0483293024 0.0695192796 0.1000000000
#>
#> $conv
#> [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
#>
#> $loss
#> [1] -0.0009277110 -0.0009277110 -0.0009277044 -0.0009276244 -0.0009272925
#> [6] -0.0009265732 -0.0009257482 -0.0009247329 -0.0009188221 -0.0009014867
#> [11] -0.0008689195 -0.0008128306 -0.0006735436 -0.0004415435 -0.0001437383
#> [16] 0.0000000000 0.0000000000 0.0000000000 0.0000000000 0.0000000000
#>
#> $fbeta
#> [1] -9.272410e-04 -9.267385e-04 -9.257101e-04 -9.236947e-04 -9.200037e-04
#> [6] -9.131254e-04 -8.995113e-04 -8.728295e-04 -8.258891e-04 -7.522001e-04
#> [11] -6.309980e-04 -4.447533e-04 -2.512276e-04 -9.964324e-05 -1.557462e-05
#> [16] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
#>
#> $nparams
#> [1] 10 10 10 10 10 10 9 7 7 7 6 5 3 3 1 0 0 0 0 0
#>
#> $betaEst
#> [1] 0.0137095845 -0.0056469817 0.0036312841 0.0063286260 0.0040426832
#> [6] -0.0010612452 0.0151152200 -0.0009465904 0.0201842371 -0.0006271410
#>
# L0
penalizedRss(bhat, R, n, penalty = "L0", lambda0 = 0.01,
lambda = c(0))
#> $beta
#> [,1]
#> [1,] 0
#> [2,] 0
#> [3,] 0
#> [4,] 0
#> [5,] 0
#> [6,] 0
#> [7,] 0
#> [8,] 0
#> [9,] 0
#> [10,] 0
#>
#> $lambda
#> [1] 0
#>
#> $conv
#> [1] 1
#>
#> $loss
#> [1] 0
#>
#> $fbeta
#> [1] 0
#>
#> $nparams
#> [1] 0
#>
#> $betaEst
#> [1] 0 0 0 0 0 0 0 0 0 0
#>