parsimony.datasets.utils

1 # -*- coding: utf-8 -*- 2 """ 3 Created on Wed Aug 7 10:35:12 2013 4 5 Copyright (c) 2013-2014, CEA/DSV/I2BM/Neurospin. All rights reserved. 6 7 @author: Edouard Duchesnay 8 @email: edouard.duchesnay@cea.fr 9 @license: BSD 3-clause. 10 """ 11 12 import numpy as np 13 from scipy import ndimage

14 #import matplotlib.pyplot as plt 15 16 17 -def corr_to_coef(v_x, v_e, cov_xe, cor):

18 """In a linear model y = bx + e. Calculate b such cor(bx + e, x) = cor. 19 20 Parameters 21 ---------- 22 v_x: Float. The variance of x, var(x). 23 24 v_e: Float. The variance of e, var(e). 25 26 cov_xe: Float. The covariance between x and e, cov(x, e). 27 28 cor: Float. The desired correlation. 29 30 Examples 31 -------- 32 >>> corr_to_coef(1, 1, 0, .5) 33 0.57735026918962573 34 """ 35 b2 = v_x ** 2 * (cor ** 2 - 1) 36 b1 = 2 * cov_xe * v_x * (cor ** 2 - 1) 37 b0 = cor ** 2 * v_x * v_e - cov_xe ** 2 38 delta = b1 ** 2 - 4 * b2 * b0 39 sol1 = (-b1 - np.sqrt(delta)) / (2 * b2) 40 sol2 = (-b1 + np.sqrt(delta)) / (2 * b2) 41 42 return np.max([sol1, sol2]) if cor >= 0 else np.min([sol1, sol2])

43

44 45 ############################################################################ 46 ## utils 47 48 49 ############################################################################ 50 ## Objects classes 51 52 -class ObjImage(object):

53 """ 54 Parameters: 55 ----------- 56 mask: ??? 57 58 std: Float. Standard deviation of latent variable 59 60 beta: Float. The coefficient of information. 61 """

62 - def __init__(self, mask=None, std=.5, beta=.5):

63 self.mask = mask 64 self.beta = beta 65 self.std = std

66

67 - def get_beta(self):

68 return self.beta

69

70 - def get_std(self):

71 return self.std

72

73 - def get_mask(self):

74 return self.mask

75 76 @staticmethod

77 - def object_model(objects, Xim):

78 """Add object variance: x_ki = coef^1/2 * o_k + (1 - coef)^1/2 * e_i 79 """ 80 sigma_o = 1 81 # Image of objects label 82 labels_im = np.zeros(Xim.shape[1:], dtype=int) 83 label = 0 84 for k in xrange(len(objects)): 85 o = objects[k] 86 label += 1 87 o.label = label 88 # A) Add object latent variable 89 mask_o = o.get_mask() 90 labels_im[mask_o] = o.label 91 obj_latent = np.random.normal(0, sigma_o, Xim.shape[0]) 92 obj_latent -= obj_latent.mean() # - 0 93 # obj_latent /= obj_latent.std() * sigma_o 94 obj_latent *= 1.0 / (obj_latent.std() * sigma_o) 95 std = o.get_std() 96 Xim[:, mask_o] = (std * obj_latent + Xim[:, mask_o].T).T 97 return Xim, labels_im

98

99 100 -class Square(ObjImage):

101 - def __init__(self, center, size, shape, **kwargs):

102 super(Square, self).__init__(**kwargs) 103 self.size = size 104 self.center = center 105 self.x_grid, self.y_grid, self.z_grid = np.ogrid[0:shape[0], 106 0:shape[1], 107 0:shape[2]]

108

109 - def get_mask(self):

110 hs = self.size / 2.0 111 mask = (np.abs(self.x_grid - self.center[0]) <= hs) & \ 112 (np.abs(self.y_grid - self.center[1]) <= hs) 113 (np.abs(self.z_grid - self.center[2]) <= hs) 114 return mask

115

116 117 -class Dot(ObjImage):

118 - def __init__(self, center, size, shape, **kwargs):

119 super(Dot, self).__init__(**kwargs) 120 self.size = size 121 self.center = center 122 self.x_grid, self.y_grid, self.z_grid = np.ogrid[0:shape[0], 123 0:shape[1], 124 0:shape[2]]

125

126 - def get_mask(self):

127 mask = np.sqrt((self.x_grid - self.center[0]) ** 2 + \ 128 (self.y_grid - self.center[1]) ** 2 + \ 129 (self.z_grid - self.center[2]) ** 2) <= self.size / 2 130 return mask

131

132 133 -class Dimaond(ObjImage):

134 - def __init__(self, center, size, shape, **kwargs):

135 super(Square, self).__init__(**kwargs) 136 self.size = size 137 self.center = center 138 self.x_grid, self.y_grid, self.z_grid = np.ogrid[0:shape[0], 139 0:shape[1], 140 0:shape[2]]

141

142 - def get_mask(self):

143 mask = np.abs(self.x_grid - self.center[0]) + \ 144 np.abs(self.y_grid - self.center[1]) + \ 145 np.abs(self.z_grid - self.center[2]) <= self.size / 2 146 return mask

147

148 149 -def get_objects_edges(objects):

150 m = objects[0].get_mask() 151 m[::] = False 152 for o in objects: 153 m += o.get_mask() 154 md = ndimage.binary_dilation(m) 155 return md - m

156

157 158 ############################################################################ 159 ## Spatial smoothing 160 -def spatial_smoothing(Xim, sigma, mu_e=None, sigma_pix=None):

161 for i in xrange(Xim.shape[0]): 162 Xim[i, :] = ndimage.gaussian_filter(Xim[i, :], 163 sigma=sigma) 164 X = Xim.reshape((Xim.shape[0], np.prod(Xim.shape[1:]))) 165 # Spatial smoothing reduced the std-dev, reset it to 1 166 if mu_e is not None: 167 X -= X.mean(axis=0) + mu_e # Also ensure null mean 168 if sigma_pix is not None: 169 X /= X.std(axis=0) * sigma_pix 170 return Xim

171

Source Code for Module parsimony.datasets.utils