parsimony.utils.utils

Source Code for Module parsimony.utils.utils

1 # -*- coding: utf-8 -*- 2 """ 3 The :mod:`parsimony.utils.utils` module includes common functions and 4 constants. 5 6 Please add anything useful or that you need throughout the whole package to 7 this module. 8 9 Created on Thu Feb 8 09:22:00 2013 10 11 Copyright (c) 2013-2014, CEA/DSV/I2BM/Neurospin. All rights reserved. 12 13 @author: Tommy Löfstedt 14 @email: lofstedt.tommy@gmail.com 15 @license: BSD 3-clause. 16 """ 17 import warnings 18 #import collections 19 from functools import wraps 20 from time import time, clock 21 22 import numpy as np 23 24 #TODO: This depends on the OS. We should try to be clever here ... 25 time_cpu = clock # UNIX-based system measures CPU time used. 26 time_wall = time # UNIX-based system measures time in seconds since the epoch. 27 time = time_cpu # TODO: Make it so that this can be changed by settings. 28 29 __all__ = ["time_cpu", "time_wall", "time", "deprecated", 30 "optimal_shrinkage", "AnonymousClass"]

31 32 #_DEBUG = True 33 34 35 -def deprecated(*replaced_by):

36 """This decorator can be used to mark functions as deprecated. 37 38 Useful when phasing out old API functions. 39 40 Parameters 41 ---------- 42 replaced_by : String. The name of the function that should be used instead. 43 """ 44 arg = True 45 if len(replaced_by) == 1 and callable(replaced_by[0]): 46 func = replaced_by[0] 47 replaced_by = None 48 arg = False 49 else: 50 replaced_by = replaced_by[0] 51 52 def outer(func): 53 @wraps(func) 54 def with_warning(*args, **kwargs): 55 string = "" 56 if replaced_by is not None: 57 string = " Use %s instead." % replaced_by 58 59 warnings.warn("Function " + str(func.__name__) + \ 60 " is deprecated." + string, 61 category=DeprecationWarning, 62 stacklevel=2) 63 return func(*args, **kwargs)

64 65 with_warning.__name__ = func.__name__ 66 with_warning.__doc__ = func.__doc__ 67 with_warning.__dict__.update(func.__dict__) 68 69 return with_warning 70 71 if not arg: 72 return outer(func) 73 else: 74 return outer 75

76 77 #@deprecated("functions.properties.Gradient.approx_grad") 78 #def approx_grad(f, x, eps=1e-4): 79 # p = x.shape[0] 80 # grad = np.zeros(x.shape) 81 # for i in xrange(p): 82 # x[i, 0] -= eps 83 # loss1 = f(x) 84 # x[i, 0] += 2.0 * eps 85 # loss2 = f(x) 86 # x[i, 0] -= eps 87 # grad[i, 0] = (loss2 - loss1) / (2.0 * eps) 88 # 89 # return grad 90 91 #def make_list(a, n, default=None): 92 # # If a list, but empty 93 # if isinstance(a, (tuple, list)) and len(a) == 0: 94 # a = None 95 # # If only one value supplied, create a list with that value 96 # if a != None: 97 # if not isinstance(a, (tuple, list)): 98 ## a = [a for i in xrange(n)] 99 # a = [a] * n 100 # else: # None or empty list supplied, create a list with the default value 101 ## a = [default for i in xrange(n)] 102 # a = [default] * n 103 # return a 104 105 #def corr(a, b): 106 # ma = np.mean(a) 107 # mb = np.mean(b) 108 # 109 # a_ = a - ma 110 # b_ = b - mb 111 # 112 # norma = np.sqrt(np.sum(a_ ** 2.0, axis=0)) 113 # normb = np.sqrt(np.sum(b_ ** 2.0, axis=0)) 114 # 115 # norma[norma < TOLERANCE] = 1.0 116 # normb[normb < TOLERANCE] = 1.0 117 # 118 # a_ /= norma 119 # b_ /= normb 120 # 121 # ip = np.dot(a_.T, b_) 122 # 123 # if ip.shape == (1, 1): 124 # return ip[0, 0] 125 # else: 126 # return ip 127 # 128 # 129 #def cov(a, b): 130 # ma = np.mean(a) 131 # mb = np.mean(b) 132 # 133 # a_ = a - ma 134 # b_ = b - mb 135 # 136 # ip = np.dot(a_.T, b_) / (a_.shape[0] - 1.0) 137 # 138 # if ip.shape == (1, 1): 139 # return ip[0, 0] 140 # else: 141 # return ip 142 # 143 # 144 #def sstot(a): 145 # a = np.asarray(a) 146 # return np.sum(a ** 2) 147 # 148 # 149 #def ssvar(a): 150 # a = np.asarray(a) 151 # return np.sum(a ** 2, axis=0) 152 # 153 # 154 #def direct(W, T=None, P=None, compare=False): 155 # if compare and T == None: 156 # raise ValueError("In order to compare you need to supply two arrays") 157 # 158 # for j in xrange(W.shape[1]): 159 # w = W[:, [j]] 160 # if compare: 161 # t = T[:, [j]] 162 # cov = np.dot(w.T, t) 163 # if P != None: 164 # p = P[:, [j]] 165 # cov2 = np.dot(w.T, p) 166 # else: 167 # cov = np.dot(w.T, np.ones(w.shape)) 168 # if cov < 0: 169 # if not compare: 170 # w *= -1 171 # if T != None: 172 # t = T[:, [j]] 173 # t *= -1 174 # T[:, j] = t.ravel() 175 # if P != None: 176 # p = P[:, [j]] 177 # p *= -1 178 # P[:, j] = p.ravel() 179 # else: 180 # t = T[:, [j]] 181 # t *= -1 182 # T[:, j] = t.ravel() 183 # 184 # W[:, j] = w.ravel() 185 # 186 # if compare and P != None and cov2 < 0: 187 # p = P[:, [j]] 188 # p *= -1 189 # P[:, j] = p.ravel() 190 # 191 # if T != None and P != None: 192 # return W, T, P 193 # elif T != None and P == None: 194 # return W, T 195 # elif T == None and P != None: 196 # return W, P 197 # else: 198 # return W 199 # 200 # 201 #def debug(*args): 202 # if _DEBUG: 203 # s = "" 204 # for a in args: 205 # s = s + str(a) 206 # print s 207 # 208 # 209 #def warning(*args): 210 # if _DEBUG: 211 # s = "" 212 # for a in args: 213 # s = s + str(a) 214 ## traceback.print_stack() 215 # print "WARNING:", s 216 217 218 -def optimal_shrinkage(X, T=None):

219 220 tau = [] 221 222 if T is None: 223 T = [T] * len(X) 224 if len(X) != len(T): 225 if T is None: 226 T = [T] * len(X) 227 else: 228 T = [T[0]] * len(X) 229 230 import sys 231 for i in xrange(len(X)): 232 Xi = X[i] 233 Ti = T[i] 234 # print "Here1" 235 # sys.stdout.flush() 236 M, N = Xi.shape 237 Si = np.cov(Xi.T) 238 # print "Here2" 239 # sys.stdout.flush() 240 if Ti is None: 241 Ti = np.diag(np.diag(Si)) 242 # print "Here3" 243 # sys.stdout.flush() 244 245 # R = _np.corrcoef(X.T) 246 Wm = Si * ((M - 1.0) / M) # 1 / N instead of 1 / N - 1 247 # print "Here4" 248 # sys.stdout.flush() 249 sum_d = np.sum((Ti - Si) ** 2.0) 250 # print "Here5" 251 # sys.stdout.flush() 252 del Si 253 del Ti 254 # print "Here6" 255 # sys.stdout.flush() 256 257 Var_sij = 0 258 for i in xrange(N): 259 for j in xrange(N): 260 wij = np.multiply(Xi[:, [i]], Xi[:, [j]]) - Wm[i, j] 261 Var_sij += np.dot(wij.T, wij) 262 Var_sij = Var_sij[0, 0] * (M / ((M - 1.0) ** 3.0)) 263 264 # diag = _np.diag(C) 265 # SS_sij = _np.sum((C - _np.diag(diag)) ** 2.0) 266 # SS_sij += _np.sum((diag - 1.0) ** 2.0) 267 268 # d = (Ti - Si) ** 2.0 269 270 # l = Var_sij / np.sum(d) 271 l = Var_sij / sum_d 272 l = max(0, min(1, l)) 273 274 tau.append(l) 275 # print "tau %f" % (l,) 276 277 return tau

278

279 280 #def delete_sparse_csr_row(mat, i): 281 # """Delete row i in-place from sparse matrix mat (CSR format). 282 # 283 # Implementation from: 284 # 285 # http://stackoverflow.com/questions/13077527/is-there-a-numpy-delete-equivalent-for-sparse-matrices 286 # """ 287 # if not isinstance(mat, scipy.sparse.csr_matrix): 288 # raise ValueError("works only for CSR format -- use .tocsr() first") 289 # n = mat.indptr[i + 1] - mat.indptr[i] 290 # if n > 0: 291 # mat.data[mat.indptr[i]:-n] = mat.data[mat.indptr[i + 1]:] 292 # mat.data = mat.data[:-n] 293 # mat.indices[mat.indptr[i]:-n] = mat.indices[mat.indptr[i + 1]:] 294 # mat.indices = mat.indices[:-n] 295 # mat.indptr[i:-1] = mat.indptr[i + 1:] 296 # mat.indptr[i:] -= n 297 # mat.indptr = mat.indptr[:-1] 298 # mat._shape = (mat._shape[0] - 1, mat._shape[1]) 299 300 301 -class AnonymousClass(object):

302 """Used to create anonymous classes. 303 304 Usage: anonymous_class = AnonymousClass(field=value, method=function) 305 """

306 - def __init__(self, **entries):

307 self.__dict__.update(entries)

308

309 - def __eq__(self, other):

310 return self.__dict__ == other.__dict__

311

312 - def __neq__(self, other):

313 return self.__dict__ != other.__dict__

314 315 316 #class EnumItem(object): 317 # def __init__(self, cls_name, name, value): 318 # self.cls_name = cls_name 319 # self.name = name 320 # self.value = value 321 # 322 # def __eq__(self, other): 323 # if not isinstance(other, self.__class__): 324 # return False 325 # 326 # return self.cls_name == other.cls_name \ 327 # and self.name == other.name \ 328 # and self.value == other.value 329 # 330 # def __hash__(self): 331 # return hash(self.cls_name) | hash(self.name) | hash(self.value) 332 # 333 # def __str__(self): 334 # return "<%s.%s: %d>" % (self.cls_name, self.name, self.value) 335 # 336 # def __repr__(self): 337 # return "EnumItem('%s', '%s', %d)" % (self.cls_name, self.name, 338 # self.value) 339 # 340 # 341 #class Enum(object): 342 # """Used to declare enumerated constants. 343 # 344 # Supposed to be similar to and used as the Enum class introduced in 345 # Python 3.4. 346 # """ 347 # def __init__(self, name, *sequential, **named): 348 ## self.__dict__["_Enum__name"] = name 349 # seq_pairs = zip(sequential, range(len(sequential))) 350 # values = {k: EnumItem(name, k, v) for k, v in seq_pairs} 351 # for k, v in named: 352 # values[k] = EnumItem(name, k, v) 353 # enums = dict(values) # , **named) 354 # for k, v in enums.items(): 355 # self.__dict__[k] = v 356 # 357 # def __setattr__(self, name, value): # Read-only. 358 # raise TypeError("Enum attributes are read-only.") 359 # 360 # def __str__(self): 361 # return "Enum: " + str(self.__dict__.keys()) 362 # 363 # 364 #Info = Enum("Info", "ok", "num_iter", "t", "f", "gap", "mu", "bound", "beta", 365 # "converged") 366 # 367 # 368 #class LimitedDict(collections.MutableMapping): 369 # """A dictionary with a constraint on the set of keys that are allowed. 370 # 371 # This class is essentially a dict, but it only allows a set of keys defined 372 # at initialisation. 373 # 374 # Parameters 375 # ---------- 376 # keys : A sequence of allowed keys. The set of keys that are allowed. 377 # """ 378 # def __init__(self, *keys): 379 # if (len(keys) == 1 and isinstance(keys[0], collections.Sequence) \ 380 # and len(keys[0]) == 0) or len(keys) == 0: 381 # self.__keys = set() 382 # 383 # elif (len(keys) == 1 and isinstance(keys[0], collections.Sequence) \ 384 # and len(keys[0]) == 1): 385 # self.__keys = set(list(keys[0])) 386 # 387 # elif len(keys) == 1 and isinstance(keys[0], collections.Iterable): 388 # 389 # self.__keys = set(keys[0]) 390 # 391 # else: 392 # self.__keys = set(keys) 393 # 394 # self.__dict = dict() 395 # 396 # def add_key(self, key): 397 # if key not in self.__keys: 398 # self.__keys.add(key) 399 # 400 # def remove_key(self, key): 401 # if key in self.__keys: 402 # self.__keys.remove(key) 403 # 404 # # Key no longer valid. Remove from dictionary if present. 405 # if key in self.__dict: 406 # del self.__dict[key] 407 # 408 # def allows(self, key): 409 # return key in self.__keys 410 # 411 # def allowed_keys(self): 412 # return list(self.__keys) 413 # 414 # def __len__(self): 415 # return len(self.__dict) 416 # 417 # def __getitem__(self, key): 418 # if key not in self.__keys: 419 # raise KeyError("'%s' is not an allowed key." % (key,)) 420 # 421 # return self.__dict[key] 422 # 423 # def __setitem__(self, key, value): 424 # if key not in self.__keys: 425 # raise KeyError("'%s' is not an allowed key." % (key,)) 426 # 427 # self.__dict[key] = value 428 # 429 # def __delitem__(self, key): 430 # if key not in self.__keys: 431 # raise KeyError("'%s' is not an allowed key." % (key,)) 432 # 433 # del self.__dict[key] 434 # 435 # def __contains__(self, key): 436 # if key not in self.__keys: 437 # raise KeyError("'%s' is not an allowed key." % (key,)) 438 # 439 # return key in self.__dict 440 # 441 # def __iter__(self): 442 # return iter(self.__dict) 443 # 444 # def clear(self): 445 # self.__dict.clear() 446 # 447 # def copy(self): 448 # info = self.__class__(list(self.__keys)) 449 # info.__dict = self.__dict.copy() 450 # 451 # return info 452 # 453 # @classmethod 454 # def fromkeys(cls, keys, value=None): 455 # info = cls(keys) 456 # info.__dict = dict.fromkeys(keys, value) 457 # 458 # return info 459 # 460 # def get(self, key, default=None): 461 # if key not in self.__keys: 462 # raise KeyError("'%s' is not an allowed key." % (key,)) 463 # 464 # return self.__dict.get(key, default) 465 ## if key in self.__dict: 466 ## return self.__dict[key] 467 ## else: 468 ## return default 469 # 470 # def items(self): 471 # return self.__dict.items() 472 # 473 # def iteritems(self): 474 # return self.__dict.iteritems() 475 # 476 # def iterkeys(self): 477 # return self.__dict.iterkeys() 478 # 479 # def itervalues(self): 480 # return self.__dict.itervalues() 481 # 482 # def keys(self): 483 # return self.__dict.keys() 484 # 485 # def pop(self, *args): 486 # if len(args) == 0: 487 # raise TypeError("pop expected at least 1 arguments, got 0") 488 # if len(args) > 2: 489 # raise TypeError("pop expected at most 2 arguments, got %d" \ 490 # % (len(args),)) 491 # 492 # if len(args) >= 1: 493 # key = args[0] 494 # default_given = False 495 # if len(args) >= 2: 496 # default = args[1] 497 # default_given = True 498 # 499 # if key not in self.__keys: 500 # raise KeyError("'%s' is not an allowed key." % (key,)) 501 # 502 # if key not in self.__dict: 503 # if default_given: 504 # return default 505 # else: 506 # raise KeyError(str(key)) 507 # else: 508 # return self.__dict[key] 509 # 510 # def popitem(self): 511 # return self.__dict.popitem() 512 # 513 # def setdefault(self, key, default=None): 514 # if key not in self.__keys: 515 # raise KeyError("'%s' is not an allowed key." % (key,)) 516 # 517 # if key in self.__dict: 518 # return self.__dict[key] 519 # else: 520 # self.__dict[key] = default 521 # return default 522 # 523 # def update(self, *args, **kwargs): 524 # info = dict() 525 # info.update(*args, **kwargs) 526 # for key in info.keys(): 527 # if key not in self.__keys: 528 # raise KeyError("'%s' is not an allowed key." % (key,)) 529 # 530 # self.__dict.update(info) 531 # 532 # def values(self): 533 # return self.__dict.values() 534 # 535 # def viewitems(self): 536 # return self.__dict.viewitems() 537 # 538 # def viewkeys(self): 539 # return self.__dict.viewkeys() 540 # 541 # def viewvalues(self): 542 # return self.__dict.viewvalues() 543 # 544 # def __format__(self, *args, **kwargs): 545 # return self.__dict.__format__(*args, **kwargs) 546 # 547 # def __eq__(self, other): 548 # if not isinstance(other, self.__class__): 549 # return False 550 # return self.__keys == other.__keys and self.__dict == other.__dict 551 # 552 # def __ge__(self, other): 553 # return self.__keys == other.__keys and self.__dict >= other.__dict 554 # 555 # def __gt__(self, other): 556 # return self.__keys == other.__keys and self.__dict > other.__dict 557 # 558 # def __hash__(self): 559 # return hash(self.__keys) | hash(self.__dict) 560 # 561 # def __le__(self, other): 562 # return self.__keys == other.__keys and self.__dict <= other.__dict 563 # 564 # def __lt__(self, other): 565 # return self.__keys == other.__keys and self.__dict < other.__dict 566 # 567 # def __ne__(self, other): 568 # keys_eq = self.__keys == other.__keys 569 # if not keys_eq: 570 # return False 571 # else: 572 # return self.__dict != other.__dict 573 # 574 # def __cmp__(self, other): 575 # keys_cmp = cmp(self.__keys, other.__keys) 576 # if keys_cmp != 0: 577 # return keys_cmp 578 # else: 579 # return cmp(self.__dict, other.__dict) 580 # 581 # def __repr__(self): 582 # return "%s(%s).update(%s)" \ 583 # % (self.__class__.__name__, 584 # self.__keys.__repr__(), 585 # self.__dict.__repr__()) 586 # 587 # def __str__(self): 588 # return "Keys: %s. Dict: %s." % (str(self.__keys), str(self.__dict)) 589