2017-12-12 19 views
0

を私はMCMCを使用してデータにODE方程式のモデルを適合する使命を帯びていますが、エラーを乗り越えることはできません。はTypeError:ロジスティック(1つの)必要な位置引数不足している:「paramsは」

TypeError         Traceback (most recent call last) 
<ipython-input-21-10d0a714b5f4> in <module>() 
    26   proposed[j] = proposed[j] + np.random.normal(0,propsigma[j]) 
    27   if (proposed[j]>0): # automatically reject moves if proposed parameter <=0 
---> 28    alpha = np.exp(logistic_loglik (proposed,time,ExRatio,sig)-logistic_loglik(par_out[i-1,],time,ExRatio,sig)) 
29    u = np.random.rand() 
30    if (u < alpha): 

<ipython-input-21-10d0a714b5f4> in logistic_loglik(params, t, data, sig) 
    3 # set up a function to return the log likelihood 
    4 def logistic_loglik(params,t,data,sig): 
--> 5  return sum(norm.logpdf(logistic(params,t),data,sig)) 
    6 
    7 # set standard deviations to be 10% of the population values 

TypeError: logistic() missing 1 required positional argument: 'params' 

問題これらの行である:これらの行の機能は、対数尤度を返すために機能を設定されている

def logistic_loglik(params,t,data,sig): 
    return sum(norm.logpdf(logistic(params,t),data,sig)) 

明らかに関与ロジスティックモデルがあります。

「params」が関数によって呼び出されていないのが私の理解ですが、これは相対的な初心者のコードを解決する方法がわかりません。

以下は、コンテキストを追加するためのコードです。

import numpy as np 
import matplotlib.pyplot as plt 
from scipy.integrate import odeint 

%matplotlib inline 

def logistic(x,t,params):  
    S, R, A = x 
    r, Nmax, delta_s, beta, gamma, delta_r, delta_a, Emax, H, MICs, MICr = params 
N = S + R 
E_s = 1 - (Emax * A**H)/(MICs**H + A**H) 
E_r = 1- (Emax * A**H)/(MICr**H + A**H) 

derivs = [r * (1 - N/Nmax) * E_s * S - delta_s * S - ((beta * S * R)/N), 
     r * (1 - gamma) * (1 - N/Nmax) * E_r * R - delta_r * R + ((beta * S * R)/N), - delta_a * A] 


return derivs 

r = 0.5 
Nmax = 10**7 
delta_s = 0.025 
beta = 10**-2 
gamma = 0.5 
delta_r = 0.025 
delta_a = 0.003 
Emax = 2 
H = 2 
MICs = 8 
MICr = 2000 

[r, Nmax, delta_s, beta, gamma, delta_r, delta_a, Emax, H, MICs, MICr] = params 
S = 9 * 10**6 
R = 10**5 
A = 5.6 
x0 = [S, R, A] 

maxt = 2000 
tstep = 1 
t = np.arange(0,maxt,tstep) 

def logistic_resid(params,t,data): 
    return logistic(params,t)-data 

logistic_out = odeint(logistic, x0, t, args=(params,)) 

time = np.array([0, 168, 336, 504, 672, 840, 1008, 1176, 1344, 1512, 1680, 1848, 2016, 2184, 2352, 2520, 2688, 2856]) 
ExRatio = np.array([2, 27, 43, 36, 39, 32, 27, 22, 13, 10, 14, 14, 4, 4, 7, 3, 3, 1]) 
ratio = 100* logistic_out[:,1]/(logistic_out[:,0]+logistic_out[:,1]) 
plt.plot(t,ratio) 
plt.plot(time,ExRatio,'h') 
xlabel('Position') 
ylabel('Pollution') 

新セル

from scipy.stats import norm 

# set up a function to return the log likelihood 
def logistic_loglik(params,t,data,sig): 
    return sum(norm.logpdf(logistic(params,t),data,sig)) 

# set standard deviations to be 10% of the population values 
sig = ExRatio/10 

# parameters for the MCMC 
reps = 50000 
nparams = 3 

# output matrix 
par_out = np.ones(shape=(reps,nparams)) 
# acceptance 
accept = np.zeros(shape=(reps,nparams)) 
# proposal standard deviations. These have been pre-optimized. 
propsigma = [0.05,20,5] 

for i in range(1,reps): 
    # make a copy of previous parameters 
    par_out[i,] = par_out[i-1,] 
    for j in range(npars): 
     proposed = np.copy(par_out[i,:]) # we need to make a copy so that rejected moves don't affect the original matrix 
    proposed[j] = proposed[j] + np.random.normal(0,propsigma[j]) 
    if (proposed[j]>0): # automatically reject moves if proposed parameter <=0 
     alpha = np.exp(logistic_loglik(proposed,time,ExRatio,sig)-logistic_loglik(par_out[i-1,],time,ExRatio,sig)) 
     u = np.random.rand() 
     if (u < alpha): 
      par_out[i,j] = proposed[j] 
      accept[i,j] = 1 

#print(sum(accept[range(101,reps),:])/(reps-100)) 


#plt.plot(par_out[:,0]) 
#plt.plot(par_out[range(101,reps),0]) 
#plt.plot(par_out[:,0],par_out[:,2]) 
plt.hist(par_out[range(101,reps),0],50) 
print('\n') 
a=np.mean(par_out[range(101,reps),0]) 

ここで、この私の最初の質問私はアドバイスしてください答えに必要なすべての情報を残しているので、もし。

すべてのヘルプは大幅にあなたの関数は次のようになります

答えて

0

をapprecieatedされるだろう:

def logistic(x,t,params) 

それは引数を取ります。しかし、あなたは2でのみ呼び出す:

def logistic_loglik(params,t,data,sig): 
    return sum(norm.logpdf(logistic(params,t),data,sig)) 

に変更し、それを:

def logistic_loglik(params,t,data,sig): 
    return sum(norm.logpdf(logistic(data, t, params), sig)) 
+0

エラーは、そのためのおかげで解決した提案された変更を実施しました。私は過去数日間、次のエラーを解決しようとしました。 'ValueError:アンパックする値が多すぎます(予想3)' = alpha.png (log.logpd(logistic(ratio、t、data、sig)); - > 5 return sum(log.pre_out [i-1、]、time、ExRatio、sig)この問題を解決することができませんでしたので、進める方法についていくつかのアドバイスを提供したいと思っていました。 –

+0

私はエラーが唯一の値を期待していると感じています。 S、R、Aの場合はこれが「x」の構造に起因する可能性があると考えているが、これを改善する方法がわからない場合や、「比率」がこれらの線のどこかに必要な場合がありますが、そこに働く –

+0

それが助けて偉大な。あなたの問題が解決されたら、回答を受け入れることができます(http://meta.stackexchange.com/questions/5234/how-does-accepting-an-answer-work)。新しい問題がある場合は、新しい質問をすることができます。 –

関連する問題