某些组之间的参数不同,而其他参数适合所有组?
Param that varies between certain groups while other params are fit across all groups?
如何允许函数的一个参数在组之间变化而其他参数适合所有组?
我正在使用 lmfit 来拟合疾病传播模型。我希望函数的指数适合所有数据点,但比例因子需要在不同组之间变化(以作为不同年份不同繁殖率的不同疾病株的代理)。
查看下面我的代码:
#### create parameters ####
params = Parameters()
params.add('tau_1', value=1.,min=0.01)
params.add('tau_2', value=1.,min=0.01)
params.add('rho_1', value=1.,min=0.01)
for i in range(0, len(sorted(variable_data.flu_season.unique()))):
season = str(sorted(variable_data.flu_season.unique())[i])
params.add('theta_' + season, value=1., min = 0.01)
#### model ####
def pop_dist_inverse_grav(params, dist, host_pop, targ_pop, flu_sea, data):
parvals = params.valuesdict()
tau_1 = parvals['tau_1']
tau_2 = parvals['tau_2']
rho_1 = parvals['rho_1']
theta_1 = parvals['theta_' + str(season)]
grav_model = theta_1 * (np.power(dist, rho_1)) / ((np.power(host_pop, tau_1)) * (np.power(targ_pop, tau_2)))
return grav_model - data
按照 M Newville 的建议,我安装了配件。代码如下:
import pandas as pd
from lmfit import minimize, Parameters, report_fit
import numpy as np
variable_data = pd.read_csv("../Data/variables_for_model_fit.csv", sep=",", header='infer')
season_list = sorted(variable_data.flu_season.unique())
#### create parameters ####
params1 = Parameters()
params1.add('tau_host', value=0.24,min=0, max = 3)
params1.add('tau_targ', value=0.14,min=0, max = 3)
params1.add('rho', value=0.29,min=0, max = 3)
# "global" parameters
for i, season in enumerate(season_list):
params1.add('theta_%d' % season_list[i], value=1000., min=1, max=1e6)
# creates theta parameters for each season
#### define model ####
def grav_dist_over_pop(dist, hist_pop, targ_pop, theta, rho, tau_host, tau_targ):
return theta**(-1) * dist**rho * host_pop**(-tau_host) * targ_pop**(-tau_targ)
#### objective function ####
def objective_1(params, dist, host_pop, targ_pop, flu_sea, data):
parvals1 = params1.valuesdict()
resid = np.zeros((len(season_data),len(season_data)))
for i, data in enumerate(data):
theta = parvals1['theta_%d' % flu_sea[i]]
rho = parvals1['rho']#_%d' % flu_sea[i]]
tau_host = parvals1['tau_host']#_%d' % flu_sea[i]]
tau_targ = parvals1['tau_targ']#_%d' % flu_sea[i]]
model = grav_dist_over_pop(dist, host_pop, targ_pop, theta, rho, tau_host, tau_targ)
resid[i, :] = model - data
return resid.flatten()
#### fit global variables ####
season_data = variable_data.sample(500)
# my dataset is huge so lmfit takes an age when fitting all of the
# theta values
host_pop = np.asarray(season_data.host_city_pop.values.tolist())
targ_pop = np.asarray(season_data.target_city_pop.values.tolist())
dist = np.asarray(season_data.distance.values.tolist())
data = np.asarray(season_data.time_to_spread.values.tolist())
flu_sea = np.asarray(season_data.flu_season.values.tolist())
result = minimize(objective_1, params1, args=(dist, host_pop, targ_pop, flu_sea, data))
report_fit(result.params)
与 lmfit
(或 scipy.optimize
和我所知道的所有其他工具)的拟合始终是 "global" 拟合,找到一组参数的最佳值用于最小化一维残差数组。需要明确的是,您 可以 优化参数以适应多个数据集(或组或季节),但您必须将问题减少到计算一维残差数组的一组参数。
对于您的问题,我建议从定义 "grav model" 开始,为单个季节或数据集的数据建模。我对这类领域一无所知(但很高兴 lmfit 可能有用!),但从你的例子来看,这可能看起来像这样(如果不是,请更正)
def grav_season(dist, host_pop, targ_pop, theta, rho, tau_host, tau_targ):
return theta * dist**rho * host_pop**(-tau_host) * targ_pop**(-tau_targ)
在我看来有 3 个自变量(dist
、host_pop
、targ_pop
)和 4 个可能是变量的参数:theta
、rho
、tau_host
、tau_targ
。同样,如果需要请更正,但对于此处的目的而言,这些细节并不那么重要。
要适合单个季节/数据集,您可以这样做
from lmfit import Parameters, minimize, report_fit
def objective(params, data, dist, host_pop, targ_pop):
pars = params.valuesdict()
model = grav_season(dist, host_pop, targ_pop, pars['theta'],
pars['rho'], pars['tau_host'], pars['tau_targ'])
return (model - data)
params = Parameters()
param.add('theta', value=1., min = 0.01)
param.add('rho', value=1., min = 0.01)
param.add('tau_host', value=1., min = 0.01)
param.add('tau_targ, value=1., min = 0.01)
result = minimize(objective, params, args=(data, dist, host_pop, targ_pop))
report_fit(result.params)
现在,对于多个季节,您可以简单地为每个季节制作一个 theta
、rho
、tau_host
、tau_targ
参数。但是,在一次拟合中一起使用所有这些数据没有多大意义。如果我理解正确,您希望 tau_host
、tau_targ
和 rho
在所有季节都具有相同的值。
为此,为全局应用的指数创建参数:
params = Parameters()
param.add('rho', value=1., min = 0.01)
param.add('tau_host', value=1., min = 0.01)
param.add('tau_targ', value=1., min = 0.01)
和每季参数:
for i, season in enumerate(sorted(variable_data.flu_season.unique())):
param.add('theta_%d' % i, value=1., min = 0.01)
param.add('rho_%d' % i, expr='rho')
param.add('tau_host_%d' % i, expr='tau_host')
param.add('tau_targ_%d' % i, expr='tau_targ')
请注意,theta_i
将独立变化,但 rho_i
等将被限制为采用变量 rho
等的值。这给出了每个参数的完整集合季节,但满足您的限制。如果看起来可能需要更多测试,这种方法允许您轻松更改其中的一个或多个以单独变化。为此,您可以说:
params['rho_7'].set(value=2.0, vary=True, expr='') # vary independently
要使用这组多季节参数,您还需要多季节数据。我不确定 dist
、host_pop
和 targ_pop
是随季节变化,还是仅 data
。我假设只有 data
随季节变化(但如果没有,应该很容易改变)。构建一个包含每个季节数据的列表,然后修改您的 objective 函数,使其看起来像:
def objective(params, season_data, dist, host_pop, targ_pop):
pars = params.valuesdict()
resid = np.zeros((len(season_data), len(season_data[0]))
for i, data in enumerate(season_data):
theta = pars['theta_%d' % i]
rho = pars['rho_%d' % i]
tau_host = pars['tau_host_%d' % i]
tau_targ = pars['tau_targ_%d' % i]
model = grav_season(dist, host_pop, targ_pop,
theta, rho, tau_host, tau_targ)
resid[i, :] = model - data
return resid.flatten()
result = minimize(objective, params, args=(season_data, dist, host_pop, targ_pop))
report_fit(result.params)
希望对您有所帮助。同样,主要建议是:
- 为您正在建模的现象创建每个季节/每个数据集的模型函数。
- 使用 'expr' 将每季参数限制为 "global values"。虽然你想要四季
具有相同的值,可能是某些参数可能随线性变化
季节,或者需要增加一些价值或其他意味着他们并不是真正独立的东西。
如何允许函数的一个参数在组之间变化而其他参数适合所有组?
我正在使用 lmfit 来拟合疾病传播模型。我希望函数的指数适合所有数据点,但比例因子需要在不同组之间变化(以作为不同年份不同繁殖率的不同疾病株的代理)。
查看下面我的代码:
#### create parameters ####
params = Parameters()
params.add('tau_1', value=1.,min=0.01)
params.add('tau_2', value=1.,min=0.01)
params.add('rho_1', value=1.,min=0.01)
for i in range(0, len(sorted(variable_data.flu_season.unique()))):
season = str(sorted(variable_data.flu_season.unique())[i])
params.add('theta_' + season, value=1., min = 0.01)
#### model ####
def pop_dist_inverse_grav(params, dist, host_pop, targ_pop, flu_sea, data):
parvals = params.valuesdict()
tau_1 = parvals['tau_1']
tau_2 = parvals['tau_2']
rho_1 = parvals['rho_1']
theta_1 = parvals['theta_' + str(season)]
grav_model = theta_1 * (np.power(dist, rho_1)) / ((np.power(host_pop, tau_1)) * (np.power(targ_pop, tau_2)))
return grav_model - data
按照 M Newville 的建议,我安装了配件。代码如下:
import pandas as pd
from lmfit import minimize, Parameters, report_fit
import numpy as np
variable_data = pd.read_csv("../Data/variables_for_model_fit.csv", sep=",", header='infer')
season_list = sorted(variable_data.flu_season.unique())
#### create parameters ####
params1 = Parameters()
params1.add('tau_host', value=0.24,min=0, max = 3)
params1.add('tau_targ', value=0.14,min=0, max = 3)
params1.add('rho', value=0.29,min=0, max = 3)
# "global" parameters
for i, season in enumerate(season_list):
params1.add('theta_%d' % season_list[i], value=1000., min=1, max=1e6)
# creates theta parameters for each season
#### define model ####
def grav_dist_over_pop(dist, hist_pop, targ_pop, theta, rho, tau_host, tau_targ):
return theta**(-1) * dist**rho * host_pop**(-tau_host) * targ_pop**(-tau_targ)
#### objective function ####
def objective_1(params, dist, host_pop, targ_pop, flu_sea, data):
parvals1 = params1.valuesdict()
resid = np.zeros((len(season_data),len(season_data)))
for i, data in enumerate(data):
theta = parvals1['theta_%d' % flu_sea[i]]
rho = parvals1['rho']#_%d' % flu_sea[i]]
tau_host = parvals1['tau_host']#_%d' % flu_sea[i]]
tau_targ = parvals1['tau_targ']#_%d' % flu_sea[i]]
model = grav_dist_over_pop(dist, host_pop, targ_pop, theta, rho, tau_host, tau_targ)
resid[i, :] = model - data
return resid.flatten()
#### fit global variables ####
season_data = variable_data.sample(500)
# my dataset is huge so lmfit takes an age when fitting all of the
# theta values
host_pop = np.asarray(season_data.host_city_pop.values.tolist())
targ_pop = np.asarray(season_data.target_city_pop.values.tolist())
dist = np.asarray(season_data.distance.values.tolist())
data = np.asarray(season_data.time_to_spread.values.tolist())
flu_sea = np.asarray(season_data.flu_season.values.tolist())
result = minimize(objective_1, params1, args=(dist, host_pop, targ_pop, flu_sea, data))
report_fit(result.params)
与 lmfit
(或 scipy.optimize
和我所知道的所有其他工具)的拟合始终是 "global" 拟合,找到一组参数的最佳值用于最小化一维残差数组。需要明确的是,您 可以 优化参数以适应多个数据集(或组或季节),但您必须将问题减少到计算一维残差数组的一组参数。
对于您的问题,我建议从定义 "grav model" 开始,为单个季节或数据集的数据建模。我对这类领域一无所知(但很高兴 lmfit 可能有用!),但从你的例子来看,这可能看起来像这样(如果不是,请更正)
def grav_season(dist, host_pop, targ_pop, theta, rho, tau_host, tau_targ):
return theta * dist**rho * host_pop**(-tau_host) * targ_pop**(-tau_targ)
在我看来有 3 个自变量(dist
、host_pop
、targ_pop
)和 4 个可能是变量的参数:theta
、rho
、tau_host
、tau_targ
。同样,如果需要请更正,但对于此处的目的而言,这些细节并不那么重要。
要适合单个季节/数据集,您可以这样做
from lmfit import Parameters, minimize, report_fit
def objective(params, data, dist, host_pop, targ_pop):
pars = params.valuesdict()
model = grav_season(dist, host_pop, targ_pop, pars['theta'],
pars['rho'], pars['tau_host'], pars['tau_targ'])
return (model - data)
params = Parameters()
param.add('theta', value=1., min = 0.01)
param.add('rho', value=1., min = 0.01)
param.add('tau_host', value=1., min = 0.01)
param.add('tau_targ, value=1., min = 0.01)
result = minimize(objective, params, args=(data, dist, host_pop, targ_pop))
report_fit(result.params)
现在,对于多个季节,您可以简单地为每个季节制作一个 theta
、rho
、tau_host
、tau_targ
参数。但是,在一次拟合中一起使用所有这些数据没有多大意义。如果我理解正确,您希望 tau_host
、tau_targ
和 rho
在所有季节都具有相同的值。
为此,为全局应用的指数创建参数:
params = Parameters()
param.add('rho', value=1., min = 0.01)
param.add('tau_host', value=1., min = 0.01)
param.add('tau_targ', value=1., min = 0.01)
和每季参数:
for i, season in enumerate(sorted(variable_data.flu_season.unique())):
param.add('theta_%d' % i, value=1., min = 0.01)
param.add('rho_%d' % i, expr='rho')
param.add('tau_host_%d' % i, expr='tau_host')
param.add('tau_targ_%d' % i, expr='tau_targ')
请注意,theta_i
将独立变化,但 rho_i
等将被限制为采用变量 rho
等的值。这给出了每个参数的完整集合季节,但满足您的限制。如果看起来可能需要更多测试,这种方法允许您轻松更改其中的一个或多个以单独变化。为此,您可以说:
params['rho_7'].set(value=2.0, vary=True, expr='') # vary independently
要使用这组多季节参数,您还需要多季节数据。我不确定 dist
、host_pop
和 targ_pop
是随季节变化,还是仅 data
。我假设只有 data
随季节变化(但如果没有,应该很容易改变)。构建一个包含每个季节数据的列表,然后修改您的 objective 函数,使其看起来像:
def objective(params, season_data, dist, host_pop, targ_pop):
pars = params.valuesdict()
resid = np.zeros((len(season_data), len(season_data[0]))
for i, data in enumerate(season_data):
theta = pars['theta_%d' % i]
rho = pars['rho_%d' % i]
tau_host = pars['tau_host_%d' % i]
tau_targ = pars['tau_targ_%d' % i]
model = grav_season(dist, host_pop, targ_pop,
theta, rho, tau_host, tau_targ)
resid[i, :] = model - data
return resid.flatten()
result = minimize(objective, params, args=(season_data, dist, host_pop, targ_pop))
report_fit(result.params)
希望对您有所帮助。同样,主要建议是:
- 为您正在建模的现象创建每个季节/每个数据集的模型函数。
- 使用 'expr' 将每季参数限制为 "global values"。虽然你想要四季 具有相同的值,可能是某些参数可能随线性变化 季节,或者需要增加一些价值或其他意味着他们并不是真正独立的东西。