/usr/share/pyshared/mvpa/clfs/glmnet.py is in python-mvpa 0.4.8-1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 | # emacs: -*- mode: python; py-indent-offset: 4; indent-tabs-mode: nil -*-
# vi: set ft=python sts=4 ts=4 sw=4 et:
### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ##
#
# See COPYING file distributed along with the PyMVPA package for the
# copyright and license terms.
#
### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ##
"""GLM-Net (GLMNET) regression classifier."""
__docformat__ = 'restructuredtext'
# system imports
import numpy as N
import mvpa.base.externals as externals
# do conditional to be able to build module reference
if externals.exists('rpy', raiseException=True) and \
externals.exists('glmnet', raiseException=True):
import rpy
rpy.r.library('glmnet')
# local imports
from mvpa.clfs.base import Classifier
from mvpa.measures.base import Sensitivity
from mvpa.misc.param import Parameter
if __debug__:
from mvpa.base import debug
def _label2indlist(labels, ulabels):
"""Convert labels to list of unique label indicies starting at 1.
"""
# allocate for the new one-of-M labels
new_labels = N.zeros(len(labels), dtype=N.int)
# loop and convert to one-of-M
for i, c in enumerate(ulabels):
new_labels[labels == c] = i+1
return [str(l) for l in new_labels.tolist()]
class _GLMNET(Classifier):
"""GLM-Net regression (GLMNET) `Classifier`.
GLM-Net is the model selection algorithm from:
Friedman, J., Hastie, T. and Tibshirani, R. (2008) Regularization
Paths for Generalized Linear Models via Coordinate
Descent. http://www-stat.stanford.edu/~hastie/Papers/glmnet.pdf
To make use of GLMNET, you must have R and RPy installed as well
as both the glmnet contributed package. You can install the R and
RPy with the following command on Debian-based machines:
sudo aptitude install python-rpy python-rpy-doc r-base-dev
You can then install the glmnet package by running R
as root and calling:
install.packages()
"""
_clf_internals = [ 'glmnet', 'linear', 'has_sensitivity',
'does_feature_selection'
]
family = Parameter('gaussian',
allowedtype='basestring',
choices=["gaussian", "multinomial"],
doc="""Response type of your labels (either 'gaussian'
for regression or 'multinomial' for classification).""")
alpha = Parameter(1.0, min=0.01, max=1.0, allowedtype='float',
doc="""The elastic net mixing parameter.
Larger values will give rise to
less L2 regularization, with alpha=1.0
as a true LASSO penalty.""")
nlambda = Parameter(100, allowedtype='int', min=1,
doc="""Maximum number of lambdas to calculate
before stopping if not converged.""")
standardize = Parameter(True, allowedtype='bool',
doc="""Whether to standardize the variables
prior to fitting.""")
thresh = Parameter(1e-4, min=1e-10, max=1.0, allowedtype='float',
doc="""Convergence threshold for coordinate descent.""")
pmax = Parameter(None, min=1, allowedtype='None or int',
doc="""Limit the maximum number of variables ever to be
nonzero.""")
maxit = Parameter(100, min=10, allowedtype='int',
doc="""Maximum number of outer-loop iterations for
'multinomial' families.""")
model_type = Parameter('covariance',
allowedtype='basestring',
choices=["covariance", "naive"],
doc="""'covariance' saves all inner-products ever
computed and can be much faster than 'naive'. The
latter can be more efficient for
nfeatures>>nsamples situations.""")
def __init__(self, **kwargs):
"""
Initialize GLM-Net.
See the help in R for further details on the parameters
"""
# init base class first
Classifier.__init__(self, **kwargs)
# pylint friendly initializations
self.__weights = None
"""The beta weights for each feature."""
self.__trained_model = None
"""The model object after training that will be used for
predictions."""
self.__trained_model_dict = None
"""The model object in dict form after training that will be
used for predictions."""
# It does not make sense to calculate a confusion matrix for a
# regression
# YOH: sorry for not clear semantics... pyvmpa is evolving,
# regressions will store RegressionStatistics within the
# confusion, so it is ok to have training_confusion
# enabled, but .regression parameter needs to be set to true,
# therefor above conditioning and tuneup of kwargs in _R
#if self.params.family == 'gaussian':
# self.states.enable('training_confusion', False)
# def __repr__(self):
# """String summary of the object
# """
# return """ENET(lm=%s, normalize=%s, intercept=%s, trace=%s, max_steps=%s, enable_states=%s)""" % \
# (self.__lm,
# self.__normalize,
# self.__intercept,
# self.__trace,
# self.__max_steps,
# str(self.states.enabled))
def _train(self, dataset):
"""Train the classifier using `data` (`Dataset`).
"""
# process the labels based on the model family
if self.params.family == 'gaussian':
# do nothing, just save the labels as a list
labels = dataset.labels.tolist()
pass
elif self.params.family == 'multinomial':
# turn lables into list of range values starting at 1
labels = _label2indlist(dataset.labels,
dataset.uniquelabels)
self.__ulabels = dataset.uniquelabels.copy()
# process the pmax
if self.params.pmax is None:
# set it to the num features
pmax = dataset.nfeatures
else:
# use the value
pmax = self.params.pmax
# train with specifying max_steps
# must not convert trained model to dict or we'll get segfault
rpy.set_default_mode(rpy.NO_CONVERSION)
self.__trained_model = rpy.r.glmnet(dataset.samples,
labels,
family=self.params.family,
alpha=self.params.alpha,
nlambda=self.params.nlambda,
standardize=self.params.standardize,
thresh=self.params.thresh,
pmax=pmax,
maxit=self.params.maxit,
type=self.params.model_type)
rpy.set_default_mode(rpy.NO_DEFAULT)
# get a dict version of the model
self.__trained_model_dict = rpy.r.as_list(self.__trained_model)
# save the lambda of last step
self.__last_lambda = self.__trained_model_dict['lambda'][-1]
# set the weights to the last step
weights = rpy.r.coef(self.__trained_model, s=self.__last_lambda)
if self.params.family == 'multinomial':
self.__weights = N.hstack([rpy.r.as_matrix(weights[str(i)])[1:]
for i in range(1,len(self.__ulabels)+1)])
elif self.params.family == 'gaussian':
self.__weights = rpy.r.as_matrix(weights)[1:]
def _predict(self, data):
"""
Predict the output for the provided data.
"""
# predict with standard method
values = rpy.r.predict(self.__trained_model,
newx=data,
type='link',
s=self.__last_lambda)
# predict with the final state (i.e., the last step)
classes = None
if self.params.family == 'multinomial':
# remove last dimension of values
values = values[:,:,0]
# get the classes too (they are 1-indexed)
rpy.set_default_mode(rpy.NO_CONVERSION)
class_ind = rpy.r.predict(self.__trained_model,
newx=data,
type='class',
s=self.__last_lambda)
rpy.set_default_mode(rpy.NO_DEFAULT)
class_ind = rpy.r.as_vector(class_ind)
# convert the strings to ints and subtract 1
class_ind = N.array([int(float(c))-1 for c in class_ind])
# convert to actual labels
classes = self.__ulabels[class_ind]
else:
# is gaussian, so just remove last dim of values
values = values[:,0]
# values need to be set anyways if values state is enabled
self.values = values
if classes is not None:
# set the values and return none
return classes
else:
# return the values as predictions
return values
def _getFeatureIds(self):
"""Return ids of the used features
"""
return N.where(N.abs(self.__weights)>0)[0]
def getSensitivityAnalyzer(self, **kwargs):
"""Returns a sensitivity analyzer for GLMNET."""
return GLMNETWeights(self, **kwargs)
weights = property(lambda self: self.__weights)
class GLMNETWeights(Sensitivity):
"""`SensitivityAnalyzer` that reports the weights GLMNET trained
on a given `Dataset`.
"""
_LEGAL_CLFS = [ _GLMNET ]
def _call(self, dataset=None):
"""Extract weights from GLMNET classifier.
GLMNET always has weights available, so nothing has to be computed here.
"""
clf = self.clf
weights = clf.weights
if __debug__:
debug('GLMNET',
"Extracting weights for GLMNET - "+
"Result: min=%f max=%f" %\
(N.min(weights), N.max(weights)))
return weights
class GLMNET_R(_GLMNET):
"""
GLM-NET Gaussian Regression Classifier.
This is the GLM-NET algorithm from
Friedman, J., Hastie, T. and Tibshirani, R. (2008) Regularization
Paths for Generalized Linear Models via Coordinate
Descent. http://www-stat.stanford.edu/~hastie/Papers/glmnet.pdf
parameterized to be a regression.
See GLMNET_C for the multinomial classifier version.
"""
_clf_internals = _GLMNET._clf_internals + ['regression']
def __init__(self, **kwargs):
"""
Initialize GLM-Net.
See the help in R for further details on the parameters
"""
# make sure they didn't specify incompatible model
regr_family = 'gaussian'
family = kwargs.pop('family', regr_family).lower()
if family != regr_family:
warning('You specified the parameter family=%s, but we '
'force this to be "%s" for regression.'
% (family, regr_family))
family = regr_family
regression = kwargs.pop('regression', None)
if regression is None:
# enforce regression by default, but regression might be used as
# a binary classifier as well, so leave it as is if it was
# explicitly specified
regression = True
# init base class first, forcing regression
_GLMNET.__init__(self, family=family, regression=regression, **kwargs)
class GLMNET_C(_GLMNET):
"""
GLM-NET Multinomial Classifier.
This is the GLM-NET algorithm from
Friedman, J., Hastie, T. and Tibshirani, R. (2008) Regularization
Paths for Generalized Linear Models via Coordinate
Descent. http://www-stat.stanford.edu/~hastie/Papers/glmnet.pdf
parameterized to be a multinomial classifier.
See GLMNET_Class for the gaussian regression version.
"""
_clf_internals = _GLMNET._clf_internals + ['multiclass', 'binary']
def __init__(self, **kwargs):
"""
Initialize GLM-Net multinomial classifier.
See the help in R for further details on the parameters
"""
# make sure they didn't specify regression
if not kwargs.pop('family', None) is None:
warning('You specified the "family" parameter, but we '
'force this to be "multinomial".')
# init base class first, forcing regression
_GLMNET.__init__(self, family='multinomial', **kwargs)
|