160
163
# lfit.calc_beta = self.calc_beta # needed for cov_beta()
166
def isestimable(C, D):
168
From an q x p contrast matrix C and an n x p design matrix D, checks
169
if the contrast C is estimable by looking at the rank of vstack([C,D]) and
170
verifying it is the same as the rank of D.
173
C.shape = (C.shape[0], 1)
174
new = np.vstack([C, D])
175
if rank(new) != rank(D):
181
Return the reciprocal of an array, setting all entries less than or
182
equal to 0 to 0. Therefore, it presumes that X should be positive in
185
x = np.maximum(np.asarray(X).astype(np.float64), 0)
186
return np.greater(x, 0.) / (x + np.less_equal(x, 0.))
188
def mad(a, c=0.6745, axis=0):
190
Median Absolute Deviation:
192
median(abs(a - median(a))) / c
197
a.shape = np.product(a.shape,axis=0)
198
m = np.median(np.fabs(a - np.median(a))) / c
204
Return the reciprocal of an array, setting all entries equal to 0
205
as 0. It does not assume that X should be positive in
208
test = np.equal(np.asarray(X), 0)
209
return np.where(test, 0, 1. / X)
213
Erase columns of zeros: can save some time in pseudoinverse.
215
colsum = np.add.reduce(matrix**2, 0)
216
val = [matrix[:,i] for i in np.flatnonzero(colsum)]
217
return np.array(np.transpose(val))
219
def rank(X, cond=1.0e-12):
221
Return the rank of a matrix X based on its generalized inverse,
225
if len(X.shape) == 2:
226
D = scipy.linalg.svdvals(X)
227
return int(np.add.reduce(np.greater(D / D.max(), cond).astype(np.int32)))
229
return int(not np.alltrue(np.equal(X, 0.)))
231
def fullrank(X, r=None):
233
Return a matrix whose column span is the same as X.
235
If the rank of X is known it can be specified as r -- no check
236
is made to ensure that this really is the rank of X.
243
V, D, U = L.svd(X, full_matrices=0)
244
order = np.argsort(D)
248
value.append(V[:,order[i]])
249
return np.asarray(np.transpose(value)).astype(np.float64)
253
A basic step function: values at the ends are handled in the simplest
254
way possible: everything to the left of x[0] is set to ival; everything
255
to the right of x[-1] is set to y[-1].
259
>>> from numpy import arange
260
>>> from models.tools import StepFunction
264
>> f = StepFunction(x, y)
268
>>> print f([[3.2,4.5],[24,-3.1]])
273
def __init__(self, x, y, ival=0., sorted=False):
278
if _x.shape != _y.shape:
279
raise ValueError, 'in StepFunction: x and y do not have the same shape'
280
if len(_x.shape) != 1:
281
raise ValueError, 'in StepFunction: x and y must be 1-dimensional'
283
self.x = np.hstack([[-np.inf], _x])
284
self.y = np.hstack([[ival], _y])
287
asort = np.argsort(self.x)
288
self.x = np.take(self.x, asort, 0)
289
self.y = np.take(self.y, asort, 0)
290
self.n = self.x.shape[0]
292
def __call__(self, time):
294
tind = np.searchsorted(self.x, time) - 1
300
Return the ECDF of an array as a step function.
302
x = np.array(values, copy=True)
304
x.shape = np.product(x.shape,axis=0)
306
y = (np.arange(n) + 1.) / n
307
return StepFunction(x, y)
309
def monotone_fn_inverter(fn, x, vectorized=True, **keywords):
311
Given a monotone function x (no checking is done to verify monotonicity)
312
and a set of x values, return an linearly interpolated approximation
313
to its inverse from its values on x.
317
y = fn(x, **keywords)
321
y.append(fn(_x, **keywords))
326
return scipy.interpolate.interp1d(y[a], x[a])
added
removed
nipy/fixes/scipy/stats/models/tools.py
