Basics

Introduction

The package implements a collection of Generalized Elastic Net (GELnet) solvers, as outlined in the following publication: https://journals.plos.org/ploscompbiol/article?id=10.1371/journal.pcbi.1004790 This vignette covers the basic usage of the gelnet package. The interface to the solvers was designed to be very flexible, by allowing the user to specify a large number of parameters. At the same time, nearly all of the parameters are given reasonable default values, making the interface easy to use if the additional flexibility is not required.

Building your first GELnet model

Let \(X\) be a samples-by-features data matrix and \(y\) be a column vector of labels. Typically, \(X\) and \(y\) are determined by the prediction task at hand, but for the purposes of this tutorial, we are going to generate them randomly:

X <- matrix( rnorm( 1000 ), 20, 50 )
y <- rnorm( 20 )

Let’s further assume that we are given a feature-feature relationship matrix \(A\). If working with genomic features, \(A\) might be the adjacency of a gene-gene interaction network. Again, we are going to generate a random matrix for the purposes of this tutorial:

A <- matrix( sample( 0:1, 50*50, repl=TRUE ), 50, 50 )
A <- A & t(A)  ## Make the matrix symmetric

We are now going to utilize the GELnet toolkit to learn a linear regression model, such that the model weights are more similar for the features that share an interaction on \(A\). As discussed in the manuscript, this can be achieved by formulating a feature-feature penalty matrix using either the graph Laplacian or \((I-D)\), where \(D\) is the graph diffusion matrix and \(I\) is the identity matrix. The gelnet package provides a function to compute the graph Laplacian from the adjacency. Here, we utilize the normalized Laplacian to keep the penalty term on the same scale as the traditional ridge regression:

library( gelnet )
L <- adj2nlapl(A)

The model can now be learned via

model <- gelnet( X, y, 0.1, 1, P = L )
## Training a linear regression model
## Running linear regression optimization with L1 = 0.100000, L2 = 1.000000
## f = 0.395022 after iteration 7

where we set the L1-norm and L2-norm penalties to 0.1 and 1, respectively. The response for new samples is computed via the dot product with the weights:

Xnew <- matrix( rnorm( 500 ), 10, 50 )
Xnew %*% model$w + model$b
##               [,1]
##  [1,] -0.218991921
##  [2,]  0.112029376
##  [3,]  0.031200391
##  [4,]  0.240312629
##  [5,]  0.176214063
##  [6,] -0.015642708
##  [7,] -0.546951047
##  [8,]  0.007920018
##  [9,]  0.138883811
## [10,] -0.328687091

Other regression problems

Linear regression is one of the three types of prediction problems supported by the package. The other two are binary logistic regression and one-class logistic regression. The latter is outlined in the following paper: http://psb.stanford.edu/psb-online/proceedings/psb16/sokolov.pdf

The package recognizes the problem type based on the class of the \(y\) argument. To train a binary predictor, we have to provide \(y\) as a two-level factor, where the first level is treated as the positive class.

y <- factor( y > 0, levels=c(TRUE,FALSE) )
model2 <- gelnet( X, y, 0.1, 1, P=L )
## Training a logistic regression model
## Treating TRUE as the positive class
## Running logistic regression optimization with L1 = 0.100000, L2 = 1.000000
## Iteration 1: f = 0.693147
## Iteration 2: f = 0.608485
## Iteration 3: f = 0.608078

If we were to score the training data using this model, we can observe that the positive samples are receiving higher scores than the negative ones

data.frame( scores= X %*% model2$w + model2$b, labels= y )
##         scores labels
## 1  -0.89901216  FALSE
## 2  -0.83224434  FALSE
## 3  -1.18637134  FALSE
## 4  -0.67708084  FALSE
## 5  -0.77444729  FALSE
## 6  -0.98702967  FALSE
## 7  -0.35902005   TRUE
## 8   0.04331193   TRUE
## 9  -0.98796815  FALSE
## 10 -0.34270976   TRUE
## 11 -0.79347006  FALSE
## 12 -1.02688260  FALSE
## 13 -0.14885098   TRUE
## 14 -0.42396270  FALSE
## 15 -0.73966627  FALSE
## 16 -0.11734756   TRUE
## 17 -0.92818971  FALSE
## 18 -0.30958543   TRUE
## 19 -0.80571106  FALSE
## 20 -0.41120076   TRUE

However, if there is class imbalance, the scores will tend to be skewed towards the class with more samples. This can be addressed by using an additional flag when training the model

model2bal <- gelnet( X, y, 0.1, 1, P=L, balanced=TRUE )
## Training a logistic regression model
## Treating TRUE as the positive class
## Running logistic regression optimization with L1 = 0.100000, L2 = 1.000000
## Iteration 1: f = 0.693147
## Iteration 2: f = 0.641172
## Iteration 3: f = 0.641149
data.frame( scores= X %*% model2bal$w + model2bal$b, labels= y )
##        scores labels
## 1  -0.4102139  FALSE
## 2  -0.3455457  FALSE
## 3  -0.7191678  FALSE
## 4  -0.1364620  FALSE
## 5  -0.2854379  FALSE
## 6  -0.5344566  FALSE
## 7   0.2123350   TRUE
## 8   0.6414655   TRUE
## 9  -0.5065721  FALSE
## 10  0.2381333   TRUE
## 11 -0.2852109  FALSE
## 12 -0.5505074  FALSE
## 13  0.4646070   TRUE
## 14  0.1366465  FALSE
## 15 -0.2343444  FALSE
## 16  0.4680499   TRUE
## 17 -0.4338957  FALSE
## 18  0.2748091   TRUE
## 19 -0.3068224  FALSE
## 20  0.1778635   TRUE

Traditionally, the loss function for logistic regression is averaged over \(n\), the number of samples. This causes every sample to make the same contribution to the loss, which is what causes the skew towards the larger class. By using the balanced flag, the problem is reformulated slightly such that the loss is averaged over the positive and negative samples separately, and then the mean of both averages is used as the overall loss.

Finally, we can build a one-class logistic regression model using just the positive samples. To train a one-class model we simply provide NULL for the \(y\) argument:

j <- which( y == TRUE )
model1 <- gelnet( X[j,], NULL, 0.1, 1, P=L )
## Training a one-class model
## Iteration 1 : f = 0.6931472 
## Iteration 2 : f = 0.5214702 
## Iteration 3 : f = 0.520937

The model can now be used as a detector that recognizes the positive samples

data.frame( scores= X %*% model2bal$w + model2bal$b, labels= y )
##        scores labels
## 1  -0.4102139  FALSE
## 2  -0.3455457  FALSE
## 3  -0.7191678  FALSE
## 4  -0.1364620  FALSE
## 5  -0.2854379  FALSE
## 6  -0.5344566  FALSE
## 7   0.2123350   TRUE
## 8   0.6414655   TRUE
## 9  -0.5065721  FALSE
## 10  0.2381333   TRUE
## 11 -0.2852109  FALSE
## 12 -0.5505074  FALSE
## 13  0.4646070   TRUE
## 14  0.1366465  FALSE
## 15 -0.2343444  FALSE
## 16  0.4680499   TRUE
## 17 -0.4338957  FALSE
## 18  0.2748091   TRUE
## 19 -0.3068224  FALSE
## 20  0.1778635   TRUE