共有パラメータと非共有パラメータの組み合わせを持つ複数のデータセットの適合方法

Question

データセット間で共有されるいくつかの変数を持つ複数のデータセットに適合しようとしています。しかし、私はこれを行うために必要なステップが不明です。以下では、私が使用しようとしているアプローチを示しました（「問題の開始はここではうまくいかない」と説明するためです）。

In this answer誰かがパラメータを共有することができます。arcoss datasetsこれは、私がいくつかの非共有パラメータを持つことができるように適合することができる方法ですか？

私はこれをどのように達成することができるか、誰かが同じ結果を達成するためのより良いアプローチを提案することができますか？ありがとう。

import numpy as np 
from scipy.stats import gamma 
import matplotlib.pyplot as plt 
import pandas as pd 
from lmfit import minimize, Minimizer, Parameters, Parameter, report_fit, Model 


# Create datasets to fit 
a = 1.99 
start = gamma.ppf(0.001, a) 
stop = gamma.ppf(.99, a) 
xvals = np.linspace(start, stop, 100) 
yvals = gamma.pdf(xvals, a) 
data_dict = {} 
for dataset in range(4): 
    name = 'dataset_' + str(dataset) 
    rand_offset = np.random.uniform(-.1, .1) 
    noise = np.random.uniform(-.05, .05,len(yvals)) + rand_offset 

    data_dict[name] = yvals + noise 
df = pd.DataFrame(data_dict) 

# Create some non-shared parameters 
non_shared_param1 = np.random.uniform(0.009, .21, 4) 
non_shared_param2 = np.random.uniform(0.01, .51, 4) 

# Create the independent variable 
ind_var = np.linspace(.001,100,100) 

# Create a model 
def model_func(time, Qi, at, vw, R, rhob_cb, al, NSP1, NSP2): 
    Dt = at * vw 
    Dl = al * vw 

    t = time 

    first_bot = 8 * np.pi * t * rhob_cb 
    sec_bot = np.sqrt((np.pi * (Dl * R) * t)) 
    exp_top = R * np.power((NSP1 - ((t * vw)/R)), 2) 
    exp_bot = 4 * Dl * t 
    exp_top2 = R * np.power(NSP2, 2) 
    exp_bot2 = 4 * Dt * t 
    return (Qi/first_bot * sec_bot) * np.exp(- (exp_top/exp_bot) - (exp_top2/exp_bot2)) 

model = Model(model_func) 

### Issues begin here ### 

all_results = {} 
index = 0 
for col in df: 
    # This block assigns the correct non-shared parameter for the particular fit 
    nsp1 = non_shared_param1[index] 
    nsp2 = non_shared_param2[index] 
    index += 1 

    params = Parameters() 
    at = 0.1 
    al = 0.15 
    vw = 10**-4 
    Dt = at * vw 
    Dl = al * vw 

    # Non-shared parameters 
    model.set_param_hint('NSP1', value = nsp1) 
    model.set_param_hint('NSP2', value = nsp2) 

    # Shared and varying parameters 
    model.set_param_hint('vw', value =10**-4, min=10**-10) 
    model.set_param_hint('at', value =0.1) 
    model.set_param_hint('al', value =0.15) 

    # Shared and fixed parameters 
    model.set_param_hint('Qi', value = 1000, vary = True) 
    model.set_param_hint('R', value = 1.7, vary = True) 
    model.set_param_hint('rhob_cb', value =2895, vary = True) 

    # One set of parameters should be returned 
    result = model.fit(df[col], time = ind_var) 

    all_results[index] = result 

Answer 1

lmfitとのフィット感は常にパラメータオブジェクトの単一のインスタンスを使用して、複数のパラメータオブジェクトを負いません。

複数のデータセットを類似モデル（おそらく同じ数学モデルであるが、各モデルごとに異なるパラメータ値を期待する）に同時にフィットさせるには、異なるコンポーネントモデルからの残差を連結する目的関数を持つ必要があります。これらのモデルのそれぞれは、Parameters（）の単一インスタンスから取得されたパラメータを持つ必要があり、各パラメータは一意の名前を持ちます。

ので、同じ機能を持つ2つのデータセットを合わせて（のパラメータ「中央」、「振幅」、および「シグマ」と、ガウスを使用してみましょう）、あなたは次に

params = Parameters() 
params.add('center_1', 5., vary=True) 
params.add('amplitude_1', 10., vary=True) 
params.add('sigma_1', 1.0, vary=True 
params.add('center_2', 8., vary=True) 
params.add('amplitude_2', 3., vary=True) 
params.add('sigma_2', 2.0, vary=True) 

などのパラメータを定義することができますあなたは可能かもしれませんが、おそらく

def residual(params, x, datasets): 
    model1 = params['amplitude_1'] * gaussian(x, params['center_1'], params['sigma_1']) 
    model2 = params['amplitude_2'] * gaussian(x, params['center_2'], params['sigma_2'] 

    resid1 = datasets[0] - model1 
    resid2 = datasets[1] - model2 
    return np.concatenate((resid1, resid2)) 

fit = lmfit.minimize(residual, params, fcn_args=(x, datasets)) 

として、第二のモデルを計算するなど、最初のデータセットと「center_2」のモデルを計算するために「center_1」、「amplitude_1」、および「sigma_1」を使用これから見ると、パラメータ値はデフォルトでは独立しています。異なるデータセットで使用するパラメータ値を共有するには、明示的に指定する必要があります（指定したリンクされた回答のように）。

あなたが同じになるようにsigma値を必要とする場合、例えば、残留機能を変更しないだろう、上記だけでパラメータの定義に：

params.add('sigma_2', expr='sigma_1') 

あなたはに追加する2つの振幅を必要とする可能性がありますいくつかの値：

params.add('amplitude_2', expr='10 - amplitude_1') 

または多分あなたは「center_2」が「center_1」よりも大きいことを確実にしたいだろうが、量によってはフィットで決定される：

params.add('center_offset', value=0.5, min=0) 
params.add('center_2', expr='center_1 + center_offset') 

これらはすべてパラメータ値を結びつける方法です。デフォルトでは、それらは独立しています。もちろん、すべてのモデルで使用されるいくつかのパラメータを持つこともできます（たとえば、パラメータ 'sigma'を呼び出して、すべてのモデルで使用します）。