Projects/laydi
Projects
/
laydi
Archived
7
0
Fork 0

Finnished the l-pls fdr

This commit is contained in:
Arnar Flatberg 2007-09-21 13:16:40 +00:00
parent 6698ebe932
commit 18f33decc7
1 changed files with 23 additions and 24 deletions

View File

@ -4,9 +4,10 @@ from scipy import zeros,zeros_like,sqrt,dot,trace,sign,round_,argmax,\
from scipy.linalg import svd,inv,norm,det,sqrtm from scipy.linalg import svd,inv,norm,det,sqrtm
from scipy.stats import mean,median from scipy.stats import mean,median
from cx_utils import mat_center from cx_utils import mat_center
from validation import pls_jkW from validation import pls_jkW, lpls_jk
from select_generators import shuffle_1d from select_generators import shuffle_1d
from engines import * from engines import pca, pls, bridge
from engines import nipals_lpls as lpls
import time import time
@ -362,8 +363,9 @@ def mahalanobis(a, loc=None, acov=None, invcov=None):
mdist = (dot(xc, covmat)*xc).sum(1) mdist = (dot(xc, covmat)*xc).sum(1)
return mdist return mdist
def lpls_qvals(a, b, c, aopt=None, alpha=.3, zx_alpha=.5, n_iter=20,center=True, def lpls_qvals(a, b, c, aopt=None, alpha=.3, zx_alpha=.5, n_iter=20,
sim_method='shuffle',p_center='med', cov_center='med',crot=True, strict=False): sim_method='shuffle',p_center='med', cov_center='med',crot=True,
strict=False, mean_ctr=[2,0,2]):
"""Returns qvals for l-pls model. """Returns qvals for l-pls model.
@ -383,32 +385,29 @@ def lpls_qvals(a, b, c, aopt=None, alpha=.3, zx_alpha=.5, n_iter=20,center=True,
""" """
m, n = a.shape m, n = a.shape
TSQ = zeros((n, n_iter), dtype='d') # (nvars x n_subsets) p, k = c.shape
n_false = zeros((n, n_iter), dtype='d') pert_tsq_x = zeros((n, n_iter), dtype='d') # (nxvars x n_subsets)
pert_tsq_z = zeros((p, n_iter), dtype='d') # (nzvars x n_subsets)
# Full model # Full model
dat = lpls(a, b, c, aopt, scale='loads') #print "Full model start"
Wcv = lpls_jk(a, b, c ,aopt, n_blocks=None, algo=algo,center=center) dat = lpls(a, b, c, aopt, scale='loads', mean_ctr=mean_ctr)
tsq_x = hotelling(Wcv, dat['W'], p_center=p_center,alpha=alpha, crot=crot, strict=strict, Wc,Lc = lpls_jk(a, b, c ,aopt)
cov_center=cov_center) #print "Full hot"
Lcv = lpls_jk(a, b, c ,aopt, n_blocks=None, algo=algo,center=center) cal_tsq_x = hotelling(Wc, dat['W'], alpha=alpha)
tsq_z = hotelling(Lcv, dat['L'], p_center=p_center,alpha=alpha, crot=crot, strict=strict, cal_tsq_z = hotelling(Lc, dat['L'], alpha=alpha)
cov_center=cov_center)
# Perturbations # Perturbations
t0 = time.time()
Vs = shuffle_1d(b, n_iter, axis=0) Vs = shuffle_1d(b, n_iter, axis=0)
for i, b_shuff in enumerate(Vs): for i, b_shuff in enumerate(Vs):
t1 = time.time() #print i
dat = pls(ac, b_shuff, aopt, 'loads', 'fast') time.sleep(.01)
Wcv = pls_jkW(a, b_shuff, aopt, n_blocks=None, algo=algo) dat = lpls(a, b_shuff,c, aopt, scale='loads', mean_ctr=mean_ctr)
TSQ[:,i] = hotelling(Wcv, dat['W'], p_center=p_center, Wi, Li = lpls_jk(a, b_shuff, c, aopt)
alpha=alpha, crot=crot, strict=strict, pert_tsq_x[:,i] = hotelling(Wi, dat['W'], alpha=alpha)
cov_center=cov_center) pert_tsq_z[:,i] = hotelling(Li, dat['L'], alpha=alpha)
print time.time() - t1
return fdr(tsq_full, TSQ, median) return fdr(cal_tsq_z, pert_tsq_z, median), fdr(cal_tsq_x, pert_tsq_x, median)