gapore 发表于 2018-12-7 10:37:06

C\C++ 提升Python效率

本帖最后由 gapore 于 2018-12-7 10:37 编辑

有没有大佬懂cython的啊?我有段程序想用cython加速一下,但是自己不懂c语言,胡乱编的效果有限

数据来源代码(Python)如下import numpy as np
import scipy.linalg as linalg

N = 2328
Y = np.random.rand(N)
L = 200   # int in range(2,N/2)
K = N - L + 1

# Embedding
X = linalg.hankel(Y,np.zeros(K))
X = X[:-K+1,:]

# Singular value decomposition (SVD)
Xt = X.T
_,sig,Vt = linalg.svd(X)
Sig = linalg.diagsvd(sig,L,L)
r = np.linalg.matrix_rank(X)
Vt1 = Vt[:r,:]
Sig1 = Sig[:r,:r]
B1 = Sig1.dot(Vt1)
B1t = B1.T

我想用最小一乘法求解矩阵 A ,使得 B1t * A = Xt ,编的python代码如下# coding = utf-8
import numpy as np
cimport cython

def FitLadReg(double[:] endog,double[:,:] exog):
   
    exog_rank = np.linalg.matrix_rank(exog)
    beta = np.ones(exog_rank)
    xstar = exog

    n_iter = 0
    diff = 10
    cycle = False
    history = []
    while n_iter < 1000 and diff > 1e-6 and not cycle:
      n_iter += 1
      beta0 = beta
      xtx = np.dot(xstar.T, exog)   
      xty = np.dot(xstar.T, endog)
      beta = np.dot(np.linalg.pinv(xtx), xty)
      resid = endog - np.dot(exog, beta)

      mask = np.abs(resid) < .000001
      resid = ((resid >= 0) * 2 - 1) * .000001
      resid = np.where(resid < 0, 0.5 * resid, 0.5 * resid)
      resid = np.abs(resid)
      xstar = exog / resid[:, np.newaxis]
      diff = np.max(np.abs(beta - beta0))
      history.append(beta)

      if (n_iter >= 300) and (n_iter % 100 == 0):
            # check for convergence circle, shouldn't happen
            for ii in range(2, 10):
                if np.all(beta == history[-ii]):
                  cycle = True
                  break

    return history

编译后调用代码如下from FitLadReg import FitLadReg
A = np.zeros((r,L))
for i in range(L):
    A[:,i] = FitLadReg(Xt[:,i],B1t)
FitLadReg函数cython编译前后效率差不多,请问要怎么改进函数代码呢?谢谢!!!

TCY 发表于 2018-12-7 12:05:57

你用numba吧……
页: [1]
查看完整版本: C\C++ 提升Python效率