2016-11-29 136 views
0

更新,以反映意見移植一個MATLAB代碼到Python

我需要端口MATLAB代碼爲蟒蛇。不幸的是,儘管在調試和谷歌搜尋方面有很多不眠之夜,但我無法讓我的代碼運行。這裏是MATLAB代碼問題:

clc; 
serie =0:50; 
serie = serie - mean(serie); 
y = cumsum(serie); 
L = length(y); 
%Calculate the detrended fluctuation short term coefficient 
npuntos = 10; 
f1=zeros(1,npuntos); 
for n1=4:16 
%Segmentation 
seg =(1:n1:L); 
%disp(length(seg)) 
yn = zeros(1,L); 
    for k = 1:length(seg)-1 
     yaux = y(seg(k):seg(k+1)-1); 
     x = 1:length(yaux); 
     A=[sum(x.^2),sum(x); sum(x),length(x)]; 
     C=[sum(x.*yaux);sum(yaux)]; 
     v=inv(A)*C; 
     a=v(1); b=v(2); 
     pen=a; 
     ord=b; 
     ytrend = pen*x + ord; 
     yn(seg(k):seg(k+1)-1) = ytrend'; 
    end 
f1(n1) = sqrt((1/seg(end)).*sum((y(1:seg(end)-1)-yn(1:seg(end)-1)).^2)); 
end 
n1=4:16; 
f1=f1(4:end); 
p1 = polyfit(log10(n1(1:end-2)),log10(f1(1:end-2)),1); 
alpha1 = p1(1); 
disp(alpha1) 

我試圖將代碼轉換成蟒蛇去如下:

import numpy as np 
 

 
data = np.arange(start=0, stop=51, step=1) 
 
data = data.transpose() 
 
data = data - np.mean(data) 
 
y = np.cumsum(data) 
 
L = len(y) 
 
# Estimate the value of alpha1 
 

 
npuntos = 12 
 
f1 = [0] * npuntos 
 
for i, n1 in enumerate(np.arange(start=4, stop=16, step=1)): 
 
    seg = np.arange(start=0, stop=L, step=n1) # Potential error 
 
    yn = [0] * L 
 
    for k in np.arange(start=0, stop=len(seg)-1, step=1): # Potential Error 
 
     yaux = y[seg[k]:seg[k + 1]-1] # Potential Error 
 
     x = np.arange(start=1, stop=len(yaux) + 1, step=1) 
 
     A = [[sum(x ** 2), sum(x)], [sum(x), len(x)]] 
 
     C = [[sum(x * yaux)], [sum(yaux)]] 
 
     v = (np.linalg.inv(A)).dot(C) 
 
     pen = v[0] 
 
     ord = v[1] 
 
     ytrend = pen * x + ord 
 
     yn[seg[k]: seg[k + 1] - 1] = ytrend 
 
    f1[i] = np.sqrt((1/seg[-1]) * sum((y[1:seg[-1] - 1] - yn[1:seg[-1] - 1]) ** 2)) 
 
n1 =np.arange(start=4, stop=16, step=1) 
 
f1 = f1[4:] 
 
xx =np.log10(n1[1: - 2]) 
 
yy=np.log10(f1[1: - 2]) 
 
print(len(xx)) 
 
print(len(yy)) 
 
p1 = np.polyfit(xx, yy, 1) 
 
alpha1 = p1[1] 
 
print(alpha1)

不幸的是,我得到類型錯誤在程序執行時,此行

p1 = np.polyfit(xx, yy, 1) 

這確實有望爲xx是長度爲9,而xx是剛剛5.使用try/catch塊如評論所說,由輸出是完全錯誤的

try: 
    f1[n1] = np.sqrt((1/seg[-1]) * sum((y[1:seg[-1] - 1] - yn[1:seg[-1] - 1]) ** 2)) 
except IndexError: 
    f1.append(np.sqrt((1/seg[-1]) * sum((y[1:seg[-1] - 1] - yn[1:seg[-1] - 1]) ** 2))) 

錯誤是固定的。

我已經通過調試器,但我不能如何錯誤。有人可以幫忙嗎? P.S-如果有人有興趣,上面的代碼應該計算出Detrended Fluctuation Analysis (DFA)

+0

這我都試過,但它給了我同樣的錯誤。你的意思是這樣的'f1 [n1] = np.sqrt((1/seg [-1])* sum((y [0:seg [-1] - 1] - yn [0:seg [-1 ] - 1])** 2))'?這也沒有用。我得到了相同的'IndexError:列表分配索引超出範圍' –

+0

您正在使用哪個版本的python? – TheBlackCat

+0

這是Python 3.5 –

回答

3

那是因爲你有npuntos = 10,使f1 = [0] * npuntos讓您f1list的大小等於10。然後你遍歷

for n1 in np.arange(start=4, stop=16, step=1): 

和訪問f1[n1]從10到15,因爲你已經使用np模塊會給你一個IndexError

UPDATE

首先,你可以使用np.zeros((5,), dtype=np.int)

其次要弄清楚IndexError的事情。就我個人而言,我不想陷入你正在解決的數學問題,所以解決方案不會是最好的。只是一個小小的變化。我相信你知道Python索引是基於零的。因爲你會開始填充你的第五元素。我不確定你想要它。你可以做enumerate(np.arange(start=4, stop=16, step=1))創建從零開始的索引到列表:

for i, n1 in enumerate(np.arange(start=4, stop=16, step=1)): 
    ... 
    f1[i] = np.sqrt((1/seg[-1]) * sum((y[1:seg[-1] - 1] - yn[1:seg[-1] - 1]) ** 2)) 

len(np.arange(start=4, stop=16, step=1))12不是你創建的大小你f110)。所以從這一點你可以創建12元素列表。

npuntos = 12 
f1 = [0] * npuntos # or np.zeros((5,), dtype=np.int) 

還是因爲它在MATLAB做,你可以做append ING值(如@nekomatic指出的),如果它是必需的。

所以你需要用

f1[n1] = np.sqrt((1/seg[-1]) * sum((y[1:seg[-1] - 1] - yn[1:seg[-1] - 1]) ** 2)) 

try/except

try: 
    f1[n1] = np.sqrt((1/seg[-1]) * sum((y[1:seg[-1] - 1] - yn[1:seg[-1] - 1]) ** 2)) 
except IndexError: 
    f1.append(np.sqrt((1/seg[-1]) * sum((y[1:seg[-1] - 1] - yn[1:seg[-1] - 1]) ** 2))) 
+3

請注意,MATLAB代碼做同樣的事情,但在MATLAB中,當您將一個索引分配給當前邊界之外的索引時,矩陣或向量會自動擴展; Python(和大多數其他語言)中的等價結構不會這樣做。 – nekomatic

+1

@nekomatic很高興知道。:)這些信息對OP來說很有用。 –

+0

@nekomatic,你如何建議我去解決這個問題?我嘗試創建一個空列表並追加值,但如果MATLAB代碼給出正確的結果,我也會得到一個錯誤 –