160412岭回归笔记

岭回归

---

import numpy as np
from scipy import linalg

x = np.random.randn(10*2).reshape(10,2)

np.savetxt("x.txt",x)
#np.loadtxt("x.txt")
x = np.hstack((np.ones((10,1)),x))
b = np.array([2,3,4])
y = np.random.randon(10)*0.5+np.dot(x,b)
np.savetxt("y.txt",y)
lam = [1,2,3,4]
xx = np.dot(x.T,x)
xy = np.dot(x.T,y)
n = 10
p = 2

tx = x[1]
ty = y[1]

txx = xx-np.outer(tx,tx)
txy = xy-ty*tx
def cvridge0(lam0):
    lamm=np.identity(p+1)*lam0
    lamm[0,0]=0
    L = np.linalg.cholesky(txx+lamm)
    lb = linalg.solve_triangular(L,txy,lower= True)
    tcoe = linalg.solve_triangular(L.T,lb)
    return (ty-np.dot(tcoe,tx))**2
np.arrary(map(cvridge0,lam))

---

def cvridge(txmty):
    txx = xx-np.outer(tx,tx)
    txy = xy-ty
    def cvridge0(lam0):
        lamm = np.identity(p+1)*lam0
        lamm[0,0] = 0
        L = np.linalg.cholesky(txx+lamm)
        lb = linalg.solve_triangular(L,txy,lower= True)
        tcoe = linalg.solve_triangular(L.T,lb)
        return (ty-np.dot(tcoe,tx))**2
    return np.array(map(cvridge0,lam))
errm = np.array(map(cvridge,x,y))
print errm

---

import numpy as np
from scipy import linalg

class Ridge(object):
    def __init__(self,x=0,y=0,inter=True,lam=0)
        np0=x.shape
        self.n = np0[0]
        self.p = np0[1]
        if inter:
            self.x = np.hstack((np.ones((self.n,1)),x))
        else:
            self.x = x
        self.y = y
        self.lam = lam
        self.lammin = 0
        self.coe = 0
        self.xx = np.dot(self.x.T,self.x)
        self.xy = np.dot(self.x.T,self.y)
        self.cverr = 0

    def ridge(self):
        lamm = np.identity(p+1)*self.lammin
        lamm[0,0] = 0
        L = np.linalg.cholesky(self.xx,lamm)
        lb = linalg.solve_triangular(L,self.xy,lower = True)
        self.coe = linalg.solve_triangular(L.T,lb)
        #return self.coe

    def cvreg(self):
        def cvridge(tx,ty):
            txx = self.xx-np.outer(tx,tx)
            txy = self.xy-ty*tx
            def cvridge0(lam0):
                lamm = np.identity(p+1)*lam0
                lamm[0,0]=0
                L = np.linalg.cholesky(txx+lamm)
                lb = linalg.solve_triangular(L,txy,lower =True)
                tcoe = linalg.solve_triangular(L.T,lb)
                return (ty-np.dot(tcoe,tx))**2
            return np.array(map(cvridge0,lam))
        errm = np.array(map(cvridge,self.x,self.y))
        self.cverr = errm.sum(axis = 0)
        self.lammin = lam[self.cverr.argmin()]
##################################################

---

x = np.loadtxt("x.txt")
y = np.loadtxt("y.txt")
lam = [1,2,3,4]
ridge1 = Ridge(x = x,y = y,lam = lam)
ridge1.cvreg()
ridge1.cverr

---

ridge1.ridge()
ridge1.coe()

---

import numpy as np
x = np.random.normal(size = (200,3))
xx = np.dot(x.T,x)
#np.identity(2)*5
print xx
xx[1:,1:] = xx[1:,1:]+np.identity(2)*5
xx

python并行

---

imoirt mkl
print mkl.get_max_threads()

$:ipcluster start -n 4 –profile=mpi

from ipyparallel import Client
#from mpi4py import MPI
import numpy as np
import time
c = Client()
view = c[:]
with view.sync_imports():
    import numpy
x = np.array(np.random.normal(size = (20000,1000)))
st = time.time()
view.scatter('a',x)
print 1,time.time()-st
view.activate()
view.block=True

st = time.time()
view.execute('z = numpy.dot(a.T,a)')
#%px z = numpy.dot(a.T,a)
print 2,time.time()-st

st = time.time()
z = np.array(view.get("z"))
print 3,time.time()-st

st = time.time()
b1 = np.sum(z,axis = 0)
print 4,time.time()-st