非线性混合模型（共9篇）
非线性混合模型 篇1
径流中长期预报一直以来是人们关注的重大课题, 在气候变化和人类活动共同的影响下, 径流的变化更加复杂化, 对其未来的精确描述十分困难。传统的径流预测方法主要有两大类:成因预测法和统计预测法[1]。成因预测法[2]多是基于研究大气环流、长期天气过程的演变规律和流域下垫面物理状况的确定性预测模型, 是径流预测研究的一个重要发展方向, 但离实用尚有较大差距。统计预测法[3]是从大量历史资料中应用数理统计的理论和方法, 寻找预报对象和预报因子之间的统计规律和关系或水文要素自身历史变化的统计规律, 建立预报模式进行预报, 常用的有时间序列法[4]、多元回归分析法[5]等。时间序列法具有原始资料搜集简便、预见期长等优点, 但它的缺点在于不能充分利用对预测量具有很大影响的气候信息和其他因素, 导致了预报的不确定性和数据的不稳定性。回归分析法虽然充分考虑了预测量对预测因子的依赖关系, 但它忽略了预测量自身的演变规律。本文将自回归和多元回归有机结合起来, 提出了一种非线性混合回归模型, 利用神经网络进行模型求解, 并应用于戴营站的年径流预报中。这种非线性混合回归模型既考虑了预测量自身演变的客观规律, 又利用了预测量对预测因子的依赖关系, 且充分考虑了水文过程的非线性关系, 模型具有一定的理论基础, 是一种值得研究的实用方法。
1 混合回归模型
1.1 多元回归模型
在影响中长期径流的众多相关因素中, 可以通过敏感性分析和弹性分析, 选择影响显著的因素, 建立多元线性回归预测模型 (Multi-Regression Model, MR) 。
设x1, x2, …, xp是p (p>1) 个线性无关的可控变量, y是表示径流量的因变量, 若它们之间的关系为:
则称为p元线性回归模型。式 (1) 中, b0, b1, …, bp, σ2都是与x1, x2, …, xp无关的未知参数, ε是随机误差。
自回归模型是在回归模型的基础上推广、导出的。
1.2 自回归模型
对于时间序列{yt}, 如果满足:
就称模型 (2) 为p阶自回归模型 (Auto Regression Model, AR) , 或者p阶自回归。式 (2) 中, a0, a1, a2, …, ap为自回归系数, ε1为模型相应的白噪声序列。
1.3 非线性混合回归模型
回归模型和自回归模型是统计模型中的两种最基本的模型, 尽管对许多随机现象, 用这两种模型都给予了很好的描述, 但对于复杂的水文现象, 单纯地使用某种模型去描述就表现出不足, 弥补单一模型不足的一个自然想法就是将两者结合起来就形成混合回归模型。文献[6]以年径流作为自回归因子, 降水、气温和用水作为回归因子, 建立了年径流预测的线性混合回归模型, 对黄河三门峡站年径流进行预测。对于理想的多输入、单输出系统, 设其输入为{x1}、{x2}、…、{x6}, 经过系统的作用后, 输出为{y}, 则线性混合回归模型为:
式中:yt为系统在t时刻的输出;a (0) , a (1) , …, a (s) 为系统的响应函数, 也称回归系数;q0, q1, …, qs为模型的阶数, 其中q0为自回归阶数, q1, …, qs为回归阶数;εt为模型的残差。
在全球气候变暖的背景下, 降水、气温等气象要素的变化引起径流的变化, 但是径流对降水、气温变化的响应并不是呈现简单的线性关系。刘昌明[7]通过大量坡面降雨试验从微观尺度上揭示了径流率随降雨强度的非线性变化现象。同时, 从宏观方面根据流域降雨径流观测资料, 采用流域水文模型的方法, 进行径流量对降雨和气温变化的响应的模拟分析, 同样证实年径流与年降水与年气温之间关系是非线性的, 影响径流变化的各因素之间是一种强耦合的非线性关系。
径流变化的影响因素与径流量之间是一个非线性映射问题, 因此, 为反应输入变量与输出变量之间的相关关系, 构造如下非线性混合回归模型:
式中:f (·) 称为非线性混合回归函数, 其中模型的阶数q0, q1, …, qs为待估计参数, q0为自回归阶数, q1, …, qs为回归阶数。
2 非线性混合回归模型的参数确定及预测效果分析
2.1 模型建立
径流形成过程是多种因素相互作用、相互联系的复杂的自然现象, 不同时间尺度下, 径流的变化受到的影响因子不一样。但总的来说, 在不考虑人类活动的情况下, 影响河流天然径流量变化的主要因素为气温、降水变化, 主要由于全球气候变化改变了与水循环密切相关的气候因子等;另外, 由于年径流的相依性, 年径流量自身也存在相关关系。因此, 流域水文系统可以看作是一个多 (因素) 输入、单 (径流) 输出系统。
在本文研究中, 输入因子为以往的年径流量{Qt-1, Qt-2, …, Qt-n}、降水量{Pt, Pt-1, …, Pt-n}和气温{Tt, Tt-1, …, Tt-n}, 输出因子为年径流量{Qt}, 由于径流变化的影响因素与径流量之间是一个非线性映射问题, 因此, 为反应输入变量与输出变量之间的相关关系, 构造如下非线性混合回归模型:
式中:f (·) 为非线性函数, q0, q1, q2是模型的阶数;εt为模型的残差。
上式 (5) 即为非线性混合回归系统水文模型的结构。
2.2 模型定阶
根据上述模型结构, 采用潮河流域戴营水文站 (1956-2000年) 的降水量、年气温和修正后的年径流来建立预测模型, 其中1956-1992年前37年的数据作为训练样本, 1993-2000年后8年作为预测检验样本。
以线性混合回归模型的阶数来确定非线性混合回归模型的阶数。线性混合回归模型的阶数确定和参数识别方法参见文献[7]。参数的估算和阶数的确定是相互联系、互为前提的, 可通过Matlab编程一次实现。计算得到自回归阶数为2, 回归阶数为2的线性混合回归预测模型如下:
因此, 影响未来径流的因子包括当年的降水、气温以及前期的降水、气温和径流量, 共有6个因子作为输入。
2.3 参数估计
根据Kolmogorov定理, 任何一个时间序列都可以看成是一个由非线性机制确定的输入输出系统的近似, 由于神经网络具有通过学习逼近任意非线性映射的能力, 故采用神经网络强大的非线性映射能力来实现式 (5) 的非线性函数f (·) 。
设p个样本集合{ (x (t) , y (t) ) |x∈Rm, y∈Rn, t=1, 2, …, p}, 式中R为映射空间。对于p个样本集合的离散空间序列, 神经网络可以完成从输入到输出的高度非线性映射, 即可以找到某种映射使得f:Rm→Rn。采用输入节点为m, 输出节点为n, 隐层节点为u的3层神经网络来实现f, 网络的输入与输出之间的关系如下:
式中:f为sigmoid函数;Φ为输入层至隐含层的激励函数;k=1, 2, …, n;t=1, 2, …, p;xi (t) 为t时刻网络的输入;
设网络总的误差小于ε, 即
根据kolmogorov定理[8], 确定非线性混合回归预测模型的BP网络隐含层节点数为13, 最终选择[1,6]网络结构进行年径流的预测。
BP网络的隐层神经元采用tansig函数作为激励函数, 输出层神经元采用logsig函数作为激励函数。采用自适应学习速率动量梯度下降反向传播traingdx函数对网络进行训练。规定的期望误差为0.01, 最多迭代次数为5000次, 学习效率为0.1。
2.4 数据处理
由于BP网络以值域为[0, 1]的S形函数作为转换函数, 因此在利用神经网络进行预测之前, 需将实际数据进行适当的变换, 规范到[0, 1]区间, 使数据主动适应神经网络, 以提高神经网络预测的准确性。文献[9]对于研究了线性变换和几种基本的非线性平滑变换, 指出对数据样本进行COS型变换的效果较好。本文采用COS型变换, 将数据样本变换到[0, 1]区间。COS型变换的基本形式为:
2.5 预测结果分析
为使训练后的网络具有较好的预测能力, 必须有足够的样本, 否则网络将无法归纳出样本集中的内在特征。但样本过多则会造成网络过度接近样本值, 从而丧失归纳和推理能力。经过多次实验, 最后选用1956-1992年共37年的数据作为训练样本, 以1993-2000年共8年的数据作为检验样本来检验模型的精度。按照上面的训练过程, 以matlab编程实现此算法, 得到1958-1992年的年径流量的拟合值如图1, 图2为戴营1958-1992年实测和模拟径流散点图。再通过网络训练得出的权值和阈值, 将后8年的检验样本输入此网络, 得到1993-2000年的各年径流预测值如图3, 图4为戴营1993-2000年实测和模拟径流散点图。
为了进一步说明非线性混合回归模型的优越性, 本文采用自回归模型 (Auto-regression Model, AR) 、多元回归模型 (Multi-Regression Model, MR) 和线性混合回归模型 (Linear Mix Regression Model, LMR) 对同样的算例进行对比计算, 4种方法的比较结果见表1。从表1所显示的结果可以看出, 基于BP网络的非线性混合回归模型的预报合格率为75%, 对径流预报效果较好, 误差变化范围缩小, 预报稳定且精度高, 充分说明了这种方法应用于径流中长期预报的有效性。
3 结 语
本文建立了以年径流量为自回归因子, 年降水和气温为回归因子的年径流预测的非线性混合回归模型, 利用BP网络进行模型求解, 并对戴营水文站的年径流进行预测, 拟合期的年径流计算值与实测值线性回归系数R2为0.825 4, 检验期的年径流计算值与实测值线性回归系数R2为0.747 3, 预测结果表明了非线性混合回归模型进行年径流预测的有效性。同时将非线性混合回归模型与自回归模型、多元回归模型和线性混合回归模型进行对比研究, 结果表明基于BP网络的非线性混合回归模型的预报合格率为75%, 对径流预报效果较好, 误差变化范围缩小, 预报稳定且精度高, 充分说明了这种方法应用于径流中长期预报的优越性。
参考文献
[1]Chow V T, Maidment D R, Mays L W.Applied Hydrology[M].New York:McGraw-Hill, 1988.
[2]Thompson S A.Simulating of climate change impacts on waterbalances in the central United States[J].Phys Geogr, 1992, 13 (1) :31-52.
[3]Freeze R A.A stochastic-conceptual analysis of rainfall runoffprocesses on a hillslope[J].Water Resources Research, 16 (2) :391-408.
[4]洪时中.非线性时间序列分析的最新进展及其在地球科学中的应川前景[J].地球科学进展, 1999, 14 (6) :559-565.
[5]王秀杰, 练继建, 费守明, 等.基于小波消噪的混沌多元回归日径流预测模型[J].系统仿真学报, 2007, 19 (15) :3 605-36-8.
[6]蒋晓辉, 刘昌明, 刘煜, 等.用混合回归模型预测黄河三门峡站年径流变化的初步研究[J].中国科学E辑, 2004, 34 (增刊) :95-102.
[7]刘昌明.黄河流域水循环演变若干问题的研究[J].水科学进展, 2004, (5) :608-614.
[8]张乃尧, 阎平凡.神经网络与模糊控制[M].北京:清华大学出版社, 2000.
[9]李存军, 邓红霞, 朱兵, 等.BP神经网络预测日径流序列的数据适应性分析与模糊控制[J].四川大学学报 (工程科学版) , 2007, 39 (2) :25-29.
非线性混合模型 篇2
形变反演模型的非线性平差
阐述了反演问题参数估计的质量要求,将反演理论中的分辨率和精度与平差中的良好统计性质作了对比分析,为使反演参数得到良好估计,提出平差处理对策.针对反演问题经常是非线性模型,提出了非线性平差的.几种方法,最后给出参数筛选的原则和相应的统计检验方法.
作 者：陶本藻 作者单位：武汉大学测绘学院,武汉市珞喻路129号,430079刊 名：武汉大学学报(信息科学版) ISTIC EI PKU英文刊名：GEOMATICS AND INFORMATION SCIENCE OF WUHAN UNIVERSITY年，卷(期)：200126(6)分类号：P207.2关键词：形变反演模型 非线性平差 系统权 假设检验
非线性混合模型 篇3
土壤地力对产量的贡献不容忽视, 唐旭等[3]得到水稻大麦轮作体系下土壤地力对大麦和水稻产量的贡献率分别为69%和75%-81%, 马常宝等[4]认为潮土区土壤地力对小麦和玉米产量的贡献分别达51.4%和54%, 国外研究也表明高地力土壤能在一定时期内维持产量的持续稳定增长[5]。因此, 对土壤地力合理评价有助于了解土地的生产力, 从而采取相应措施提升土壤地力以尽可能发挥土壤的增产潜能。土壤地力评价比较容易的是通过施肥和不施肥产量对比或依据地力调查来评估[6,7], 目前较多的是采用更为复杂的模型方法进行评价[8,9,10], 骆东奇等[11]对土壤肥力评价指标和方法进行了总结, 但基于剖面上土壤理化性质的多元混合线性回归法则不多见, 该法考虑到不同层次上理化性质间的差异性和相互作用, 引入随机效应使得预测将更加精准[12]。
世界上最早的水稻种植区发现于湖南澧县[13], 而在道县则发现了迄今为止最早的古栽培稻谷[14], 因此, 本研究在两地分别采集了高中低三种产量的稻田土壤剖面样品, 分层对其理化性质进行测定, 然后建立稻田地力与理化性质间的多元线性回归模型和多元混合线性模型并对其进行比较分析, 同时比较不同土壤深度层次划分采样时所建立的多元混合线性模型优劣, 用于更好指导取样分析, 从而期望为稻田地力评估和指导农业生产实践提供理论参考。
1 材料与方法
1.1 研究区概况
本文选取的稻田样点分别位于湖南澧县和道县。澧县属亚热带湿润季风气候, 平均海拔59 m。该县春温秋和, 夏热冬寒, 雨量充足, 年均气温16.5℃, 年降水量1213 mm, 日照年平均1770 h, 无霜期265 d。全境以水稻土, 红壤和潮土为主, 地貌特征主要为平原和岗地, 农业种植以水稻、棉花、油料作物为主。道县也属亚热带湿润季风气候区, 平均海拔181 m。境内热量丰富, 无霜期长, 雨量充沛, 光照充足, 光能潜力大。境内年均气温18.6℃, 无霜期302 d, 年降雨量1523 mm, 年日照时数1514 h, 年蒸发量1586 mm。境内季风气候明显, 冬季盛吹北风, 气温较低;夏季盛吹偏南风, 气温较高;春秋两季温度适中, 因此四季分明。
1.2 样品采集与测定
利用土壤柱状圆筒采样器 (Eijkelkamp, 荷兰, 内径9 cm) 和剖面挖掘方法于2011年11月分别采集了两县稻田的高中低产50 cm深土壤剖面样品各一个。澧县高产田为稻-油轮作, 水稻品种为徽两优3号, 油菜品种为沛油520, 一季稻高中低产田稻谷产量分别为10.65、7.5、5.25 t/hm2。道县稻田种植两季, 高中低产田单季产量分别为7.5、6.0、5.25 t/hm2。
对于每一个土壤剖面, 按从上到下每隔5 cm取原状土环刀 (100 cm3) 并相应取该层次土壤剩余样品, 环刀样用于测定土壤两个物理性质, 即土壤饱和导水率 (Ks) 和容重 (BD) , 对应层次的剩余土样风干后, 分别过2 mm和0.25 mm筛用于测定土壤理化性质, 包括物理性质:机械组成;化学性质:阳离子交换量 (CEC) 及单个的钾 (K+) 、钙 (Ca2+) 、钠 (Na+) 和镁 (Mg2+) 离子浓度, p H, 电导率 (EC) ;养分特征:总氮 (TSN) , 总磷 (TSP) , 总钾 (TSK) , 有机碳 (SOC) , 速效磷 (AVP) , 速效钾 (AVK) 和碱解氮 (AVN) 。
CEC采用氯化钡提取-快速滴定法[15], AVN采用碱解扩散法[16], SOC利用碳氮分析仪测定, Ks采用双环刀恒定水头法, 土壤理化性质测定主要参照刘光崧[17]主编的《土壤理化分析与剖面描述》, 机械组成利用比重计法测定, 采用美国制, 即砂粒 (SAND) :0.25-2.0 mm;粉粒 (SILT) :0.25-0.002 mm;粘粒 (CLAY) <0.002 mm。
1.3 数据处理与统计
1.3.1 数据处理
本文依据所选稻田水稻产量情况, 将稻田地力定量为不同数值, 即高产=3;中产=2;低产=1。若认为土壤理化性质在同一土壤剖面中不同深度层次上相互独立, 将所有理化性质与稻田地力作相关性分析, 选取显著的影响因子, 利用逐步回归法构建多元线性回归模型。
若考虑在同一剖面中, 不同深度层次同一土壤理化性质间具有关联性, 则对土壤深度分成四组进行考虑, 分别为:10层, 即0-50 cm按5 cm进行分层;5层:即0-50 cm按10 cm进行分层;3层:即为0-20 cm, 20-30 cm, 30-50 cm;分2层, 即按犁底层划分为0-20 cm和20-50 cm。由于分层情形不同, 每层中所含的实测样品数则不同, 若分10层, 则每层含有6个样品数据;分5层, 则每层含有12个数据;分3层, 则对应层样品数分别为24, 12, 24;分2层, 则对应层样品数分别为24和36。利用每层中样本数据进行相关性分析, 选出在每层上都具有显著性的土壤理化性质, 然后分别建立四个多元混合线性回归模型并比较, 从而得到预测效果最好的分层采样方法用于指导取样分析。
1.3.2 统计方法
依据不同假设, 本研究分别建立多元线性回归模型和Laird和Ware[12]提出的多元混合线性模型来预测稻田地力并比较。
多元线性回归模型一般表达式:
式中:Y为因变量向量, X为自变量向量, β为回归系数向量, ε为残差向量。
多元混合线性模型一般表达式:
式中:i为分组, Yi表示第i组因变量向量, Xi为第i组固定效应向量, βi为第i组固定效应系数向量, Zi为第i组随机效应向量, λi为第i组随机效应系数向量, εi为第i组残差向量。
多元混合线性模型, 其实质是对多元线性回归模型的残差项作进一步的拆解, 进而引出随机效应, 它能够获得更准确的预测效果。多元混合线性模型常用于对聚类或非独立的分层数据以及同一个体不同时间段的重复测量数据进行拟合[18], 该模型同时纳入了对不同聚类或分层数据产生的固定效应和随机效应影响, 固定效应都是相同的, 不同只是随机效应。随机效应主要源于从总体中取样时带来的波动性变化, 对样品数据的影响具有非系统性和难以预测性, 而固定效应则相反, 对样品数据具有系统性和可预见性特征[19]。
利用均方根误差 (RMSE) 和调整的确定系数 (R2adj) 来评估多元线性回归模型和多元混合线性模型预测优劣, RMSE值越小, 表明模型预测越准确。R2adj位于0到1之间, 其值越接近1, 模型拟合越准。而针对不同分组构建的不同多元混合线性回归模型比较则采用方差分析, 其评判依据是赤池信息准则 (Akaike information criterion, AIC) , 贝叶斯信息准则 (Bayesian information criterion, BIC) 和残差 (ε) , 这些值越小, 表明模型拟合越好。评判指标的计算公式分别表示为:
式中:ym和yp分别为样品观测值和预测值, n为观测样品数。R2为模型确定系数, M为变量个数。k为模型中未知参数个数, L为模型最大似然函数值,
利用R统计软件 (http://www.r-project.org/) 进行相关性分析、作图和模型构建。多元混合线性模型中选用的软件包为lmer 4.0, 建模前进行相关性分析并分析不同分组情形下相关系数分布, 然后提取显著的固定效应 (P<0.05) 和随机效应进行建模, 不同分组情况下查表得到的对应相关系数临界值分别为:10层:0.755 (n=6) ;5层:0.553 (n=12) ;3层:0.553 (n=12) ;2层:0.396 (n=24) , 即认为相关系数大于该临界值时表明相关性显著, 多元混合线性模型构建过程中会自动计算出相应的AIC、BIC和ε。
考虑到Ks具有较强空间变异性, 通常符合对数正态分布[20], 因此本文对Ks进行了常用对数转换, 即K10=log10 (Ks) 。
2 结果与分析
2.1 土壤性质
高中低产稻田土壤剖面样品所有理化性质特征见表1, 土样质地以粉壤土和壤土为主, 分别占45%和32%。在三种不同地力稻田土壤中, K10都表现出较强的空间变异性, 且在高产田中变异系数最大, 但K10的平均值在高产田中最小, AVP在高产田中也表现为强变异 (CV>100) , 而其他性质在不同地力土壤中都属于中等变异 (10<CV≤100) 。高产田中的SILT, CLAY, SOC, AVK, p H, CEC, Ca2+和K+较中、低产田都相对高, 而BD, TSP, TSK, AVN, AVP, EC, Mg和BCSR在中产田中平均值较高, 只有SAND在低产田中最大, TSN含量在高中低产田中差异不明显。
注:K10:饱和导水率;BD:容重;CLAY/SILT/SAND:粘粒/粉粒/砂粒含量;SOC:有机碳;TSN:总氮;TSP/TSK:总磷/总钾;AVN/AVP/AVK:碱解氮/速效磷/速效钾;p H:p H值;Ca2+/K+/Mg2+/Na+/CEC:钙/钾/镁/钠/阳离子交换量;EC:电导率;BCSR:盐基饱和度;CV:变异系数 (%) 。
2.2 因子筛选
从多元线性回归模型中的土壤理化性质与稻田地力间相关系数可知 (表2) , 影响地力的显著土壤性质中, Ca2+和CEC影响最大, 相关系数都为0.59 (P<0.01) , 其次为K+和AVK, 相关系数分别为0.44 (P<0.01) 和0.39 (P<0.01) , 最后为p H和CLAY, 其他性质都不显著。Na+, EC, BCSR, K10和SAND表现出与地力的负相关性, 其他性质都为正相关性。
注:YIELD_N:稻田地力;K10:饱和导水率 (cm/d) ;BD:容重 (g/cm3) ;CLAY/SILT/SAND:粘粒/粉粒/砂粒含量 (%) ;SOC:有机碳 (%) ;TSN:总氮 (%) ;TSP/TSK:总磷/总钾 (g/kg) ;AVN/AVP/AVK:碱解氮/速效磷/速效钾 (mg/kg) ;p H:p H值;Ca2+/K+/Mg2+/Na+/CEC:钙/钾/镁/钠/阳离子交换量 (cmol/kg) ;EC:电导率 (u S/cm) ;BCSR:盐基饱和度 (-) 。*:P<0.05;**:P<0.01。
从多元混合线性模型中的土壤理化性质与稻田地力间相关系数随土壤深度变化可知 (图1) , 化学性质方面, Ca2+, CEC和p H表现为从上层至下层的同步增长趋势且都显著, 并且Ca2+和CEC在20-30cm处呈现出于地力最强的相关性, 而K+则呈现同步下降趋势且都显著, Mg2+, EC和BCSR在不同分层上对稻田地力影响都不显著, Na+在土壤深度10-30 cm处表现出与稻田地力较强的负相关性。养分特征中只有AVP, AVK和TSP在不同层次上与地力关系都紧密 (P<0.05) , 从土壤剖面上层到下层, AVK都呈现正相关, 分布比较均匀, 对于TSP则在上层相关性比在下层强, 而AVP在20-30 cm处与地力相关显著性不明显。AVN、TSN和TSK与地力的相关系数在各深度层次上都较低且在不同层次上的相关性呈正负交替变化, SOC在20-30 cm处呈现较强的正相关性, 再往下则变为负相关。物理性质中, K10在上层与稻田地力为显著正相关, 20 cm以下则呈现显著负相关。BD与K10有着较为相似的相关系数分布趋势。CLAY在下层对稻田地力影响较大且为正相关, 而SAND则在下层表现出与地力较强的负相关。通过对图1的分析, 相关性从上层至下层转变剧烈且对稻田地力影响较大的性质有K10, SOC, Na+和BD, 这些性质在多元线性回归建模中都没有表现出较好的相关性, 但是经过分层处理后则在不同层次上会有显著性变化出现, 表明这些性质受分层影响, 建立混合效应模型时将作为随机效应, 而其他显著变量如Ca2+, CEC, p H, K+, AVP, AVK, TSP, CLAY和SAND则将作为固定效应进行模型构建。
2.3 模型分析
建立的多元线性回归模型如下:
该模型只选用了四个变量, 分别为p H、K+、Ca2+和TSP, 能解释稻田地力总变异的73%, 这说明还有一些未知的因素未考虑到而影响预测精度, 方程中K+和TSP回归系数较高, 表明模型对其变化比较敏感。模型中未选用显著变量CEC和AVK, 而选用p H和TSP, 主要是由于Ca2+和CEC, 以及K+和AVK之间的显著相关性, 在多元逐步回归计算中进行了筛选, 同时考虑到其他较为重要的性质。
不同分组情形下建立的多元混合线性模型对稻田地力的预测比较得出, 土壤深度分为3层, 即将土壤剖面按0-20、20-30、30-50 cm分层采样所构建的多元混合线性模型预测效果最好, 此时AIC和BIC值都最小 (表1) , R2adj达到0.85, 表示该模型能解释稻田地力的85%。与多元线性回归比较时, 多元混合线性模型的RMSE较小且R2adj较高, 表明多元混合线性模型预测比多元线性回归预测更准确。
土壤分3层构建的稻田地力多元混合线性模型包括固定效应 (表4) 和随机效应 (表5) , 选取的固定效应包括AVP、AVK、TSP、Ca2+、K+、CEC、p H、CLAY和SAND;随机效应包括K10、SOC、Na+和BD。总变异中由Na+和BD所带来的变异比值较大, 分别占模型总变异值的64.5%和9.7%, 表明这两个变量对稻田地力的多元混合线性模型会产生更大的随机影响。
注:Ca2+/K+/Mg2+/Na+/CEC:钙/钾/镁/钠/阳离子交换量 (cmol/kg) ;p H:p H值 (-) ;EC:电导率 (u S/cm) ;BCSR:盐基饱和度 (-) ;AVN/AVP/AVK:碱解氮/速效磷/速效钾 (mg/kg) ;TSN:总氮 (%) ;TSP/TSK:总磷/总钾 (g/kg) ;SOC:有机碳 (%) ;K10:饱和导水率 (cm/d) ;BD:容重 (g/cm3) ;CLAY/SILT/SAND:粘粒/粉粒/砂粒含量 (%) ;虚线为对应相关性显著临界值:即10层:0.755 (n=6) ;5层:0.553 (n=12) ;3层:0.553 (n=12) ;2层:0.396 (n=24) 。
利用所建立的多元混合线性模型和多元线性回归模型观测值与预测值作箱线图比较分析 (图2) , 多元混合线性模型对稻田地力预测值相对多元线性回归模型都较为集中, 尤其是对中、低地力的预测更好, 而在高地力预测方面, 简单模型预测值虽然表现较为集中, 但存在较多异常值, 进而影响了模型整体估测精度。
3 讨论
3.1 土壤理化性质对稻田地力的影响
土壤地力表征着土壤的基础潜能, 反映了土地的出产力。土壤肥力是土壤地力的重要属性, 高肥力土壤意味着高地力。在高地力土壤上往往在相对较少的人为管理或施肥下都能获得较高作物产量[3]。本研究将高中低产稻田的地力情况定量化, 即高产田代表高地力土壤, 中产田代表中等地力, 低产田代表低地力土壤。许多研究者认为土壤自身理化性质及肥力特性对作物产量起重要作用[21,22], 本研究发现高产田中土壤平均SILT, CLAY, SOC, AVK, p H, CEC和Ca2+含量较高, 而K10最低。SILT和CLAY代表土壤颗粒大小和含量, 细质地土壤有助于保持土壤养分进而提升土壤肥力。从图1中看出, K10在表层 (0-20 cm) 表现为与地力较强的正相关性, 底层为负显著相关性。表层K10越大, 土壤导水性越强, 说明土壤空隙较充分且土壤较为疏松, 这有助于养分运输和分配, 而底层若K10越大, 则土壤将难以保持充足养分和水分, 表明高地力土壤应在表层有较高的导水性而在底层则导水率应小。
注:Int:截距;AVP/AVK:速效磷/速效钾 (mg/kg) ;TSP:总磷 (g/kg) ;Ca2+/K+/CEC:钙/钾/阳离子交换量 (cmol/kg) ;p H:p H值 (-) ;CLAY/SAND:粘粒/砂粒含量 (%) 。
注:Int:截距;K10:饱和导水率 (cm/d) ;SOC:有机碳 (%) ;Na+:钠离子 (cmol/kg) ;BD:容重 (g/cm3) 。
图2多元线性回归模型和多元混合线性回归模型对稻田地力的预测值与观测值的比较Fig.2 Comparison between the observed and predicted soil productivity by the multiple linear regression model (left panel) and the mixed-effects linear regression model (right panel)
养分性质方面, SOC表征着土壤有机质, 吴槐泓等[23]通过定位试验认为有机质能指示土壤肥力水平并与水稻产量成正相关。本文从多元线性回归模型中的相关性分析发现SOC与地力的相关性很弱 (表2) , 但经过分层处理后, 发现在土壤20-30 cm深度处SOC与地力表现较强的正相关性。SOC在土壤中的含量一般随深度而逐渐降低, 但本研究认为犁底层处较高的SOC含量有助于提升稻田地力。AVK在高产田中含量高说明提高土壤中有效钾含量有助于提升地力 (图1) , 蒋德安等[24]认为适当供钾可以显著提高水稻谷粒产量, 这表明土壤中钾素的不足将限制着土壤地力, 但从图1发现, 土壤TSK并不与地力呈现显著相关, 只有土壤中AVK和K+增加, 才能显著提升稻田地力。
化学性质方面, 土壤深度10 cm以下, p H与地力都呈现显著正相关 (图1) , 表明底层高p H有助于增强稻田地力, 杨曾平[25]通过室内盆栽试验认为以石灰岩风化物发育的水稻土肥力比河流沉积物和第四纪红色粘土发育的水稻土高。CEC和Ca2+含量在高产田中表现比较突出, 与稻田地力呈显著正相关, 表明高地力稻田对阳离子的需求量比较大, 尤其是钙离子。本研究发现低产田中砂粒含量较大, 这跟砂粒不利于养分和水分的保持有关。通过对比高中低产稻田土壤理化性质分析发现, 稻田地力不受土壤总氮的影响, 这可能与目前的施肥制度有关, 通过施肥提高产量已成为目前农田中最常用的手段, 因此不同水稻产量稻田中总氮量已趋于饱和, 可能实际上对于高地力土壤本不需要施入当前用量的氮肥, 多余的氮则将导致氮肥浪费和环境污染。研究还发现某些土壤性质在考虑分层影响后, 与地力的相关性在某些层次上呈现显著相关性, 如Na+在20-30cm处表现较强的负相关性, K10, BD和SOC表现为表层较强的正相关而底层较强的负相关性 (图1) , 表明这些土壤性质对地力的影响受分层影响。
3.2 稻田地力模型
多元线性回归得到的稻田地力预测模型只包含了p H, K+和Ca2+以及TSP四个变量, 变量相对较少, 但模型预测精度只有73%, 并且多元线性回归模型要求样本数据需具有独立性, 然而考虑到土壤剖面样品在不同层次间可能存在相互影响, 如SOC随深度减小变化, 上层SOC与下层SOC可能具有较好的相关性, 因此建立了预测稻田地力的多元混合线性模型, 即同时选取固定效应和随机效应, 固定效应选取与地力具有显著相关性但基本不随土壤分层而发生正负相关性变化的土壤性质, 而随机效应则选那些因分层而导致的显著性及相关性正负变化的性质, 随机效应在不同层次上对地力总体预测表现并不一致, 由于混合效应考虑到了更多变异来源, 对重复数据或聚类或分层数据具有较好的应用前景[26,27], 理论上多元混合线性模型将提高预测精度。本文得到的多元混合线性模型将预测精度从多元线性回归模型的73%提高到了85%, 但多元混合线性模型选用了比较多的固定效应, 这将增加实验分析难度。因此, 应对于某些变量进一步筛选, 如TSP, K+, AVP, SAND和CLAY对地力的影响在不同层次影响其实并不都完全显著, 同时, 可能还有一些更好的变量未考虑, 如王远敏[28]通过田间控制试验认为稻田施硅肥有利于产量增加。
4 结论
本研究发现不同地力稻田土壤总氮差异不明显, 低产田中砂粒平均含量偏高而其他土壤理化性质都偏低。通过混合效应模型中相关性分析筛选出13个预测稻田地力的因素, 即固定效应9个, 包括速效钾、速效磷、总磷、钙离子、钾离子、阳离子交换量、p H值、粘粒和砂粒含量, 其中阳离子交换量与钙离子浓度对稻田地力影响最显著且正相关。随机效应4个, 即饱和导水率, 有机碳, 钠离子和容重, 其中钠离子浓度产生的随机效应最强。建立的多元线性回归模型对钾离子浓度和总磷含量变化比较敏感, 本文建立的稻田地力混合效应模型较多元线性回归模型预测性能从73%提高到了85%。
人口预测模型的非线性动力学研究 篇4
人口问题不仅是20世纪我国所面临的最重大的问题之一，而且在新世纪中将继续存在。无论是对我国目前经济发展状况的认识，还是对未来经济发展的预测，人口问题的研究都具有十分重要的意义。
众所周知，人口的增长并不是按比例线性增长的，也就是说人口问题是非线性动力学问题。在非线性系统中，可能会出现分岔现象，分岔理论为我们研究人口问题提供了一种新的方法。分岔和混沌是非线性系统特有的现象，而实际上我们所处的社会经济系统中的绝大部分都是非线性，因此可以运用分岔和混沌理论来研究。本文将分岔和混沌理论用以研究人口问题，提出了人口的跨临界分岔。
一、模型的建立
设在时刻t(以年为单位)时的人口数量为p(t)，则人口的增长率为，再假设出生率b(p,t)和死亡率d(p,t)不随时间变化，则由马尔萨斯人口理论可得：
附图
其中：d[,1]=b(p,t)-d(p,t)>0，一般而言出生率大于死亡率。
显然马尔萨斯模型存在着重大的缺陷，它没有考虑到物种之间的竞争、自然界的平衡和人文环境因素，即生活资料及空间的局限、人与人的竞争、生产力水平、文化水平以及传统观念等因素。
生活资料及空间的局限和人与人的竞争会导致冲突，对人口的增长起到制约的作用。统计规律显示，从人口p(t)中随机抽取一个人，他与其他人冲突的概率与人口总数成正比，即：冲突次数=k[,1]p(t)，那么个人冲突的总次数为。而制约作用也是随着总次数的增长而增长的，同总次数成正比，以δ[,1](p,t)表示制约作用，则
附图
可以认为影响作用与生产力发展水平α[,2]成反比，与文化发展水平α[,3]成反比，与中国传统观念α[,4]成正比，与人口数量p(t)成正比，用δ[,2](p,t)表示这种影响，则：
附图
假设在未来的若干年内，各地区的人口比例保持不变，则对于某个地区来说，机械迁徙的人口数量与该地区经济发达程度α[,5]成正比，与迁出区的人口数量成正比，由于可以近似地认为各地区的人口比例不变，所以可以认为与该地区的人口数量成正比，则：
附图
二、模型分析
附图
运用稳定性分析方法，可以很容易地判断其稳定性。
附图
图中，实线表示稳定，虚线表示不稳定。由分岔图可以看出该模型发生跨临界分岔。由分岔图可以看出，在α的某个领域内，当α>0时，人口数量是稳定增长的;而当α<0时，从α系数的物理含义上看，说明死亡率大于其他影响因素系数之和，但是此时由于先前α>0，人口数量已经有所增加，所以能保持人口数量的稳定，同时我们也可以看出只有短期内保持α<0，人口数量才会稳定。
将再分岔出周期8解、16解等等，大约在3.56994时，进入混沌区，如图2所示。
三、模型的求解
设初始时刻t[,0]的人口为p[,0]，将(5)式分离变量后两边取定积分，进行求解可得：
附图
只要求出α和β，就可以运用(9)式作人口预测。理论上，只要将两年的年份与人口数量，即：t[,1]、P(t[,1])、t[,2]、P(t[,2])代入(9)式，则可得一个二元方程组，通过求解此二元方程组，便可以得到α和β值，但是此二元方程组为超越方程组，求解十分困难。
生态学家通过大量的统计，认为不鼓励也不限制生育时的α值为0.029，将α代入(9)式则可求得：，可以看出β会随着所代入的t的不同而变化，不再是一个常数，而是一个变量;另一方面，我国实行计划生育，限制了出生率，并且模型(5)中考虑了生活资料和空间的局限性、人与人的竞争、生产力水平、文化水平以及中国传统观念等因素的影响，因此α的值不再是0.029。
基于上述原因，本文采用数值解法，即：确定适当的.α和β，使得：
附图
其中：p[,i]是第i年的人口数量。则
附图
根据(11)式求解十分复杂，可以用Matlab或c语言编写一段程序，在一定范围内搜索α与β值，使得(10)式成立。
四、实证分析
本文就河北省1952年～人口数量进行实证分析和检验。众所周知，在建国后近50年中，我国经历了3年自然灾害，80年代开始实行了计划生育政策，这些都会改变动力系统的特征，因此不能将这49年的数据看成是同一动力系统的，而应该将它们分成不同的动力系统：1952年～1959年、1966年～1979年、1980年～1989年、1990年～20。3年自然灾害造成了人口数量的骤减，不仅改变了动力系统的特征，而且还造成了1960年～1965年数据的不稳定，应予以剔除;80年代后实行计划生育，但是1989年的**引起了大量人口的机械迁移，因此90年代视为另一个动力系统。
表1数据计算和检验(单位：万人)
附图
由表1可以看出，使用本文所提出的方法，能够较为精确地预测出人口的数量。四个阶段的α均大于0.029，这表明生活资料和空间的局限性、人与人的竞争、生产力水平、文化水平、中国传统观念、人口机械迁移以及实行计划生育等因素的综合影响，使人口趋于增多。相对人口增长率，即，总体上越来越小，80年代的攀升主要是由于人口迁移的影响，从整体上看人口自然增长率(一般小于人口增长率)也是逐年降低的。
【责任编辑】顾岚
【参考文献】
1方亚玲.《对人口模型的研究》.《山西煤炭管理干部学院学报》，年第2期.
2罗警山.《人口模型中的分岔研究》.《系统工程与电子技术》，1990年第10期.
非线性混合模型 篇5
一、原理与方法
线性混合效应模型的一般形式表示为
其中 分别是已知和随机设计矩阵。
β是未知的固定因素向量, 对所有个体参数都适用。
是未知的随机因素向量, 不同的个体可以有不同的随机向量。
是随机误差.假设随机效应服从标准正态分布。即:
是组间方差-协方差矩阵, 是组内协方差矩阵。所谓的重复测量可以是指在同一时间点上用不同的测量方法或者工具测量多次, 也可以是在不同的时间点, 用同一种测量方法或者工具测量一次, 也可以是前两情况的结合。假设在m个时间点上的重复测量情况, 那么这些测量来自于以下正态分布协方差矩阵:
或者可以将之简写为 m×m阶组内协方差矩阵)
这种情况, 因为是用同一测量方法或者工具, 所以不同时间点之间的测量是相关的。
线性混合模型通过调整随机设计矩阵 和定义协方差矩阵 来调试模型。常见运用中, 对均衡完整资料的重复测量设计通常假定 的边际协方差或者是非结构化的 (Unstructured, 即任意两点间的相关不等或不完全相等) 或者是等相关结构化的 (Compound Symmetry, 也称复合对称性, 即任意两时点的相关是相等的) 。后者可以直接通过指定 而得到。
注:1代表当满足右边括号里的条件时, 左边的值存在 (乘以1)
表示第i个和第j个协变量之间的欧几里德 (Euclidean) 距离, 由输入数据的观测值决定
在线性混合效应模型中, 通过指定G的一般协方差结构和任一个R的常用结构, 可以得到整体的方差-协方差结构。表1的第一和第二种是常见的协方差结构, 同时, 可以看到后面三种协方差结构因为考虑了距离等因素, 更有空间概念。这里笔者将引入本文的关注点, 实际实验中, 如果是非等间距的实验, 那么方差-协方差结构应该要考虑空间距离, 而不是简单地选择非结构化的或等相关结构化的协方差结构。下面, 笔者通过实例就这点做些说明。这些协方差结构的选择都可以通过SAS中的PROC MIXED中的选项实现。
二、实例分析
例:欲了解动物实验中某药物用药计量对牛体重的影响。实验入组40头牛 (16头母牛, 16头公牛) , 完全随机区组设计。试验中, 各受试的牛将在基线, 第15天, 第21天, 第29天和第43天分别测试体重。描述统计分析见表2。
考虑到基线体重也是体重变化的一个重要因素, 我们将引进基线体重作为线性混合效应模型中的协变量, 同时考虑固定因素:组别, 性别和测试时间, 随机因素:实验编号和入组号, 可参考以下SAS程序:
简单说明:random语句放的是上面提到的随机因素。Repeat语句放的是重复测量的变量基于每个测试动物个体的每个测试时间 (day) , 程序中&type.处可以填写模型所需的协方差结构类型, UN, CS, SP (POW) (Day) 等等。SAS输出结果中‘Fit Statistics’显示模型拟合AIC值, 值越小拟合结果越好。分别用不同的协方差结构进行测试, 得到下表3:
通过比较, 我们发现常规使用的UN和CS的协方差结构已经不是最优选择, 而空间幂函数
(SP (POW) ) 的AIC值最小, 从而被作为协方差结构新的选择, 模型调整为:
参考文献
[1]SAS for Mixed Model
非线性混合模型 篇6
对显式时间积分算法如中心差分法,当质量和阻尼矩阵为对角矩阵时,得到的动力平衡方程为非耦合,因而可容易地实现并行化.然而,显式算法为有条件稳定难以有效地应用于结构动力分析[3],隐式算法为无条件稳定而受到青睐,但不易直接实现并行.由于非线性分析每一次迭代需要求解方程组,隐式时间积分的最耗时部分为方程组求解.由于有效刚度矩阵为对称和正定,可采用并行线性求解技术实现并行计算[4].
结构动力分析的并性算法研究一直很活跃,Noor在流水线向量机实现了显式积分法的并行计算[5],Chiang等[6]在共享储存MIMD环境下对隐式Newmark类并行算法进行了研究,Noor等[7]在多处理器计算机中设计了基于混合格式的有限元算法.以往的并行研究主要在向量机和共享储存MIMD并行计算机上开展.近年来,基于网络机群环境的分布式MIMD并行计算系统成为并行计算机发展的主流方向之一[8].本文在网络机群并行计算环境下,综合隐式和显式时间积分技术,对结构非线性动力反应分析进行混合时间积分算法的并行研究,算法用可移植的MPI软件开发环境实现,通过算例验证了该算法的正确性和有效性.
1 Newmark时间步算法
结构动力分析的控制方程可表示为[9]
式中,M为质量矩阵,C为阻尼矩阵,K为刚度矩阵,a,v和d分别为加速度,速度和位移向量,f(t)为荷载向量.求解这个问题需对时间进行积分,可通过对时间的离散化实现.求解式(1)的最广泛使用的直接时间积分方法是Newmark类方法[3]该类方法是基于第(n+1)时间步的平衡,即在时间(tn+1=tn+Δt)的平衡来表示
式(2)中的位移和速度向量可写成以下有限差分形式
参数β和γ确定Newmark类算法,当和γ=可得到Newmark恒定加速度法.该法为无条件稳定和二阶精度,是结构动力分析的有效方法.
用前一步计算的已知值dn,vn和an,来确定式(2)～式(4)的3个未知量dn+1,vn+1和an+1.由式(3)和式(4)求an+1和vn+1,并代入式(2),dn+1可由下式求得
式中,为有效刚度矩阵Keff.除非计算步大小改变,在线性分析的所有计算步中,有效刚度矩阵为常数.在非线性分析的每一个时间步中,有效刚度矩阵和线性算法需要修改.非线性分析有效刚度矩阵可表示为,式中KT为切线刚度矩阵.有效残余力向量r可由下式得到
式中p为内力向量.可用Newton-Raphson类的迭代法求解这些非线性方程.
对显式时间积分,常采用中心差分法.vn和an的差分格式可表示为
将该差分格式代入动力平衡方程,可得到
式中,p为内力向量.若M和C为对角矩阵,求解dn+1不需要分解.由于中心差分法不需要组集整体矩阵,并可以在自由度级进行处理,这样可直接实现并行.与隐式算法相比较,每个时间步的处理器间通信要求很少[3].
2 混合时间积分
混合时间积分算法是在隐式和显式网格划分中分别采用隐式和显式算法.最初用于串行求解结构与介质的相互作用[10].采用区域分解技术,该算法可有效地开发出并行性.区域分解时,将有限元网格划分为与处理器数相同的子区域数.每个子区域分配给一个处理器,相应子区域的计算在各处理器中进行.这样在隐式时间积分算法中求解的整体方程组可划分为较小的子集并在各处理器中独立求解.子区域的界面采用并行显式算法进行求解.
把每个子区域中的节点分为内部节点(内部单元节点)和界面节点(位于子区域间界面上的节点),采用隐式-显式算法,界面节点采用并行显式算法进行显式积分.当界面值为已知,利用得到的已知值限制界面节点自由度,内部节点可隐式积分.由于在界面采用显式积分,利用隐式算法积分的各子区域被完全分开.这样,隐式算法就能用于局部子区域.
为简便说明这种扩展的混合时间积分算法的信息流,采用一维情况进行描述,如图1所示.水平轴表示一个一维有限元网格,垂直轴表示时间.节点间的信息流用带箭头直线表示.界面节点信息流用虚线表示.在tn+1时更新节点i,需要节点i-1,i,i+1在前一时间tn时的信息.图中实线表示子区域内的信息流,其节点采用隐式积分.
为开发混合时间积分算法的并发性,而不影响隐式算法的无条件稳定,对显式和隐式算法采取不同的时间步.隐式算法的时间步通常取为显式时间步的整数倍.但这种整数倍时间步实现不是直接的,对一维网格如图2所示,隐式积分时间步Δti取为显式积分时间步Δte的2倍.从图中的显式算法计算流可以看出,一些内部节点需要在每个显式时间步更新.在内部节点更新信息的唯一方法是采取显式算法对这些节点进行积分.如图2,若Δte为Δti的一半,每个子区域的一个内部节点需要进行积分,大部分实际情况中,Δti为Δte的10～20倍.当取Δti=10Δte,则在t+Ate时9个内部节点需要积分,在t+2Δte时8个内部节点需要积分,如此类推.这样影响算法的并行实现.本文采取弱耦合的方法,即隐式和显式节点完全分开积分.对N个时间步界面节点采用显式算法积分,隐式时间步Δti取为显式时间步Δte的N倍.然后,利用界面节点作为边界节点,在各个子区域采用隐式算法进行积分.这样就可避免在隐式和显式区域涉及多个时间步实现的问题.基于这种弱格式的算法如下.
将有限元网格划分为子区域
(1)确定内部和边界节点
(2)基于稳定准则确定显式积分时间步大小,即确定Δte
(3)建立隐式积分时间步大小,即Δti=N×Δte(N为整数)
(4)设置时间步大小Ate
(5)从1到N循环
更新时间
采用并行显式算法积分界面节点
结束循环
(6)约束界面节点
(7)设置时间步大小Δti
(8)各处理器中每个子区域进行隐式时间积分
(9)如(时间步数
3 混合时间积分算法的并行实现
在基于区域分解的并行混合时间积分算法中,将位于界面的节点划分为主边界节点和从属边界节点.每个界面节点在一个子区域中为主边界节点,而在共享该节点的所有其他子区域中为从属边界节点.非界面节点称为内部节点.不包含边界节点的单元为内部单元,而包含边界节点的单元为边界单元[11].并行混合时间积分算法的数据划分如图3所
■主边界节点□从属边界节点〇内部节点
将有限元网格划分为一定数目的子区域,每个子区域的主边界节点和从属边界节点即可确定,并行混合时间积分首先是利用显式时间积分方法对边界节点积分.在每个处理器中对各自主边界节点进行显式时间积分,将计算出这些节点的位移与相邻处理器进行通信.包含主边界节点的单元内力向量的计算也在其各自处理器中进行,并将部分值发送给其相邻的对应部分,同时从相邻处理器接收相应的部分内力向量.这样就完成了有限元网格的边界节点的时间积分,而不会出现负载不平衡.这些界面节点的时间积分将继续到N个显式时间步.
下一步是在各处理器中计算局部子区域的动力反应,利用显式积分得到的数值约束边界节点,在局部子区域中采用Newmark隐式算法进行时间积分来实现动力反应的计算.这种实现,可避免整体方程组的求解.方程组的求解只局限于子区域中.这样就显著减小通信时间和计算时间,并且局部子区域的矩阵带宽比整体矩阵要小得多.
4 数值算例
本文研制的非线性动力分析混合时间积分并行算法的数值计算是在DELL工作站机群上进行.该机群由4台双CPU的DELL工作站通过100.0Mbps以太网连接而成,共有8个CPU (2.4 GHz Xeon chips,512KB cache),每个节点内存1.0 GB.每台工作站都有真实的IP地址,使用消息传递接口(MPI)[12].操作系统为Red Hat Linux9.0.
如图4所示为1/4的圆筒形壳,L=1000mm,两直边固定.受垂直均布荷载p=2.0kN/m2,施加方式如图4.壳体厚度2 mm,E=2.06×105 MPa,泊松比v=0.3,屈服应力fy=235MPa,质量密度为7.8×103kg/m3.考虑几何非线性.
将壳划分为12×24结构网格,采用8节点的壳单元.通过采用时间步为0.2×10-3 s,来验证算法的稳定性和精度.用隐式积分算法(Newmark算法)和本文混合时间积分算法计算节点A位移时间历程反应,计算结果如图5所示.可以看出本文并行混合时间积分算法的结果与Newmark算法的结果非常吻合,表明本文并行混合时间积分算法是有效的,即具有很好的稳定性和精度.
为测试算法的性能,加密有限元网格以增加问题大小(自由度).为保证解的精度,采用严格的收敛准则.并行算法的加速比如图6和图7所示,本文算法与Newmark算法的CPU时间(自由度为4710)如表1所示.可以看出,对相同规模的问题本文算法比隐式算法快,相应的加速比高.并且算法的计算性能随着问题大小的增加而提高,显示出本文并行混合时间积分算法具有良好的性能.
5 结论
综合隐式和显式时间积分技术,对混合时间积分的并行算法进行了研究.该算法采用区域分解技术.在并行混合时间积分算法中,通过在子区域界面采用显式时间积分技术,引入并发性.这就易于分解子区域并能在每个处理器中独立求解.由于不需要并行求解器,使求解方程工作显著减小.方程求解是在局部子区域中进行,整体矩阵的半带宽显著减小.数值算例表明本文算法具有良好的性能.
参考文献
[1] Belytschco T,Yen H J,Mullen R.Mixed methods for time integration.Comput Meth Appl Mech Engng,1979,17/18: 259-275
[2] Hajjar JF,Abel JF.Parallel processing for transient nonlinear structural dynamics of three dimensional framed structures using domain decomposition.Comput & Struct, 1988,30:1237-1254
[3] Sziveri J,Topping BHV.Transient dynamic nonlinear analysis using MIMD computer architectures.J Comput Civil Engng,2000,14(2):79-91
[4] Bennighof JK,Wu JY.An unconditionally stable parallel transient algorithm with multi-step independent subdomain computation.Comput Syst Engng,1991,2:217-230
[5] Noor AK.Lambiotte.Finite element dynamics analysis on CDs Star 100 computer.Comput & Struct,1979,10(1): 7-19
[6] Chiang KN,Fulton RE.Structural dynamics methods for concurrent processing computers.Comput & Struct,1990, 36:1031-1037
[7] Noor AK,Peter JA.Partitioning strategy for efficient nonlinear finite element dynamic analysis on multiprocessor computers.Comput & Struct,1989,31:795-810
[8] Danielson KT,Namburu RR.Nonlinear dynamic finite element analysis on parallel computers using Fortran 90 and MPI.Adv Engng Software,1998,29:179-86
[9]付朝江,张武.非线性动力有限元重叠区域分裂的隐式并行算法.计算力学学报,2006,23(6):783-788(Fu Chaojiang, Zhang Wu.An implicit parallel algorithm for nonlinear dynamic finite element analysis with overlapped domain decomposition. Chinese Journal of Computational Mechanics, 2006,23(6):783-788(in Chinese))
[10] Rama Mohan Rao A.A time integration algorithm for nonlinear dynamic analysis.Adv Engng Software,2005,33: 261-271
[11] Day DM,Bhardwaj MK.Mechanism free domain decomposition. Comput Methods Appl Mech Engrg,2007,192: 763-776
非线性混合模型 篇7
在工程计算中常会遇到非线性代数方程组的求解问题。解非线性方程组F(X)=0,通常方法有两大类:一类属于线性化方法,即用一线性代数方程组求近似逼近非线性代数方程组,由此构造一组递推公式,用于逐次逼近所求的根。这类方法主要有牛顿迭代法及各种改进,用Newton方法求解,初始点的选取适当否决定了算法的收敛性及收敛速度;另一类方法是把方程组所求解问题化为多元函数的极小值的等效问题求解。这类方法主要度法及其各种改进。
梯度法的优点是程序简单,每次迭代所需的计算量小,所要求的存储量变少,而从一个不太好的初始近似出发,还是能收敛于局部极小值点。这个方法最大的最大缺点是收敛速度慢。(线性速度),尤其在解的附近,常为了提高一点精度而需要付出很大的代价。牛顿法具有收敛快和自校正等优点,其收敛阶大于等于2,但初始值选取不当易导致求解失败。结合两算法优点,先采用梯度法搜索出X*,当F值降到一定程度以后,再改用修正Newton法进行若干步搜索。这样,把梯度法与牛顿法相结合,就能使两个方法相互取长补短,在迭代初始值不太好的情况下也能保证收敛性,加快收敛速度。
2 问题描述
设有如下非线性方程组的求解问题:
x∈D,D奂Rn有界,F(f1(X),f2(X),…fn(X))T=0,其中,fi:Rn→R(i=1,2,…n)是连续可微的实值函数。若存在X*∈D使F(X*)=0,则X*称为方程组(1)的解。
用最速下降法解方程组(1)时,通常将方程组(1)构造一目标函数:
于是非线性方程组(1)的解就是上式(2)的零极小值点,反之亦然。
3 混合算法
步骤1:给定初值X0、最速下降法允许误差E、控制常数C、修正Newton法计算精度e1、e2及简化步m值。
步骤2:用最速下降法搜索出X*。
步骤3:用修正Newton法再进行若干步搜索,若未找出比X*更好的点,结束。否则,将修正Newton法找到的更好的点替代X*,返回步骤1。
4 实例与程序设计过程
验证方法的可行性与有效性,举例进行数值实验。
陈文青:一种求解非线性方程组的混合优化算法及程序实现
4 结论
本程序采用C语言编写完成,对求解非线性方程组,读者只需要根据具体所求解方程组适当修改即可,本程序的编写中,还可采用更好的方式实现,优化代码。从计算结果可以看出,把梯度法与牛顿法相结合,能使两个方法相互取长补短。先采用最速下降法搜索出X*,当F值降到一定程度以后,再用修正Newton法进行若干步搜索,若未找出比X*更好的点,结束。否则,将修正Newton法找到的更好的点替代X*,这样在迭代初始值不太好的情况下也能保证收敛性,同时又加快收敛速度,特别当初始值X0选取离解X*比较远时,效果显著。
参考文献
[1]李庆杨.非线性方程组的数值解法[M].北京:科学出版社,1999.
[2]Ortega J M,Rheinboldt W G.Iterative Solution of Nonlinear Equations in Several Variables[M].New York:Academic Press,1970.
[3]郑咸义.数值计算方法与FORTRAN语言[M].北京:电子工业出版社,1986.
[4]角仕云.实用科学与工程计算方法[M].北京:科学出版社,2000.
非线性混合模型 篇8
《非线性编辑》课程是教育技术学专业的一门专业必修课, 是一门应用性、操作性、实践性较强的课程。本课程重点学习Adobe Premiere非线性编辑软件, 主要讲授非线性编辑的基本知识、视频编辑技术、特效的设置与调节、字幕及音乐的合成与输出, 与Cool3d、Photo Shop、After Effects、Particle Illusion 3.0等软件的综合使用, 培养学生影视节目编辑、节目包装制作的能力。其主要任务是及时让学生进行选题、撰写文稿、分镜头、拍摄、采集、编辑、音乐、音效、解说、特技、合成等一系列的实践操作环节, 即能够独立使用非线性编辑系统完成节目制作, 也能够通过网络与他人共享素材协同工作。
九江学院作为一所地方本科院校, 在教育技术学、广播电视编导、影视多媒体等专业开设《非线性编辑》课程, 重在培养学生的软件应用能力、创新精神、团队协作及问题解决的能力。本人作为课程的主讲教师之一, 积极尝试在教学实践中应用不同的教学模式来指导教学, 并寻找一种能适合本校学生学习的模式。经过几年的教学实践探索, 发现混合式学习在本门课程中的效果是最优的, 下面笔者就结合实际的经验谈谈混合式学习在《非线性编辑》课程教学中的应用。
二、混合式学习的原理
混合式学习就是“一种将面授教学与基于技术媒介的教学相结合而构成的学习环境”[1], 它借助面授教学与基于技术媒介的教学这两种学习模式的优势来重新组织和实施教学活动, 以达到提高教学效率的目标。麦森和莱恩尼[2] (Mason&Rennie) 则进一步扩展了格林汉姆的定义, 认为“混合式学习是技术、场所、教学方法的多方面融合”, 而不仅仅是教学组织形式的结合。在国内教育技术界, 北师大何克抗教授认为, 混合式学习“就是要把传统学习方式的优势和E-learning的优势结合起来”。
新媒体联盟2015地平线报告高等教育版中, 在短期趋势中也提出了混合学习的广泛应用。“混合式学习” (Blended Learning) 是在线学习和面授相结合的学习方式, 从本质上来讲混合学习是把传统学习方式的优势和网络化学习 (E-learning) 的优势结合起来进行教学。也就是说, 在教学中既要发挥教师的引导、启发、监控教学过程的主导作用, 又要充分体现学生作为学习主体的主动性、积极性和创造性。它借助面授与网络这两种学习模式的优势来重新组织教学资源、实施学习活动, 把这两者结合起来, 使两者的优势互补, 以达到最佳的学习效果。
三、混合式学习在《非线性编辑》课程中的应用现状
Moodle是基于建构主义教育理论而开发的课程管理系统, 是一个免费的开放源代码的软件[3]。《非线性编辑》课程就是基于Moodle平台建立的, 根据课程教学的需求, 可以添加教学日历、点到、论坛、测验、资源、投票、问卷调查、作业、聊天、专题讨论、Wiki、互动评价等模块, 界面部分截图如图1所示, 并在此基础上进行二次开发增加新的页面补充教学资源, 以便学生更好地利用Moodle平台进行混合式学习, 教师将上课所用的教学文件、课件、教学素材、拓展资料、学生毕业设计作品等上传到Moodle平台, 学生可以用自己的帐号登录, 既可以提前预习课程内容也可以在课后复习与拓展学习, 遇到问题可以在线与教师、同学交流, 解决疑问。同时, 为每个学生提供1G空间, 将制作的作品上传到个人博客, 教师观看并给出成绩。这些在《非线性编辑》课程教学的过程中, 都发挥了积极的作用。
由于网络技术、移动终端及网络教育的迅猛发展, 学习者的学习环境正在发生着翻天覆地的变化。基于移动互联网的信息技术在教育领域中的应用越来越广。混合式学习的基本形式是师生面对面与学生在线学习的混合, 而网络平台或者APP恰恰突破了时空的限制, 打破了传统课堂讲授的单一模式。在《非线性编辑》课程具体的教学中还没有充分体现出混合式学习的特点, 没有充分利用好网络的共享性、开放性、交互性好特点, 资源建设也很缺乏。基于以上的问题, 至此, 对混合式学习在《非线性编辑》课程中的应用有了新的想法和做法。
四、混合式学习在《非线性编辑》课程中的优化措施
与传统教学不同的是, 在混合式教学中, 既强调教师在教学中的指导性, 同时也注重发挥学生学习的主动性。尤其突出的是, 强调师生在整个学习过程中的相互沟通、交流和互动, 这就需要采取与传统教学不同的教学策略、设计方案和课件工具, 通过便于移动学习的短小精悍的微课及具有展示性、互动性和形象性的交互式视频课件, 将会对《非线性编辑》课程的教学起到一定的促进作用。
1. 制作微课。
随着我校的信息化迅速发展, 学生个人拥有的信息终端 (如手机、平板电脑、上网本、笔记本电脑等) 逐步普及, 学生利用碎片化时间及非正式学习方式的增多, 按以往的将一节课45分钟的课堂教学录像视频, 放到教学资源平台上供学生学习, 已无法顺应移动互联网时代学习者的“注意力模式”, 难以满足学生的实际需求, 造成了教育资源建设与应用过程中的极大浪费。因此, 以阐释某一知识点为目标, 以短小精悍的在线视频为表现形式, 以学习或教学应用为目的的在线教学视频的微课诞生了。微课的时长一般在10分钟以内, 这也正好是学生注意力集中的最佳时间。由于《非线性编辑》课程中讲解软件应用的课占到了总课时的一半, 对于转场效果、特效应用及插件的使用等具体的知识点就可以用专业摄像机摄制、屏幕录制的方式将其录制成微课, 然后上传到课程网站, 这样的视频既短小精悍又占用空间小, 学生既可以在网站进行重复播放学习也可以下载到移动终端进行学习, 有利于学生随时随地个性化的学习, 图2、图3是用专业摄像机在录播教室录制的微课, 而图4则是通过屏幕录制的方式录制的微课。
2. 交互式视频课件的制作与应用。
《非线性编辑》课程要让学生充分实现混合式学习的模式, 除了面授所需的资源外还需要大量的网络资源, 微课为学生提供了可视化的直观的视频资料, 为使《非线性编辑》的教学更加生动, 吸引学生主动学习, 正在尝试制作交互式视频课件。交互式视频课件 (Interactive Video Courseware, 简称IVC) 是一种基于多种软件技术和工具开发出来的数字化、模块化和结构化的多用途电子资源包[4]。从结构上看, 它是以视频摄制与编辑为核心, 以互联网为发布平台, 结合各种互动课件工具组合而成的一种综合性、多样化的大型课件形式[5]。其主要功能是利用各种基于计算机和互联网的数字化技术 (Digital Technology) 来呈现、展示和传播特定内容的信息或内容, 能够最大程度地提高资讯传播的速度、范围和效率。交互式视频课件可分为交互设计类、语音动画类、视频编辑类和整合发布类四大类。交互设计类的课件可用Raptivity、Articulate Studio等软件制作生成Flash格式, 从而便于跨平台使用和传播;而语音动画可以用Character Builder和Inter Phonic来实现, 图5是就是运用Character Builder及Inter Phonic来实现课件中模拟教师语音授课和语音合成, 为学生创建了一个尽量接近现实教学环境的虚拟情景。本课件的主要内容是《非线性编辑》课程介绍, 课件中人物为Susan (body) 模型, 讲解的声音是通过Inter Phonic将文字用语音转换为Wave格式的诵读语音文件。课件是用PPT制作好后导入的, 根据课程内容用Character Builder可快速定制和生成Flash格式的动画人物, 并对人物的表情、运作进行个性化定制开发。对于介绍性比较强的知识点可以利用Character Builder及Inter Phonic软件来制作介绍性及互动性的课件。在以后的教学中, 还将引进虚拟视频展示、远程视频面对面等资源扩充混合式学习。
3. Moodle《非线性编辑》课程平台的优化。
目前基于Moodle的《非线性编辑》课程, 直接以Moodle自带的功能在运行, 界面的美观、布局的合理及功能的实用性方面比较欠缺, 在不断的优化过程中改进现有功能, 积极引入网上答题、现场直播等新技术。
五、混合式学习在《非线性编辑》课程中的展望
根据学习金字塔[6]的实验结果, 不同教学方法下, 学生学习24小时之后对新学习内容平均保持的百分比如图6所示。我们的观点是对学生多感官刺激、多媒体观看以及主动参与, 课堂不只是教师一人讲授, 学生能向其他同学讲授他所学习的知识, 解释这些知识, 小组讨论并能当堂运用这些知识, 或者通过移动终端在线或直播观看, 都会大大提高学生对所学信息的保持效果。这种混合式教学模式如图7所示。
在操作性、实践性很强的《非线性编辑》课程教学中实施混合式学习, 希望学生在以下三个方面有所突破。
1. 学生主动学习。
教学设计合理, 课件制作精良, 时间在10分钟以内, 能引导学生主动学习, 调动学生的学习积极性, 促进学生学习能力提高, 教学信息清晰明确的微课;设计上突出强调教学过程的完整性, 如教学内容展示、师生互动、作业、测验和评估等的交互式视频课件。这些优质教育资源上网, 学生通过网络可自主学到最优质的课程, 掌握学习的主动权。
2. 实现翻转课堂。
教师对教学业务熟练、课程内容讲解富有逻辑性、思维清晰、条理分明, 这是学生比较喜欢的, 但是随着微课、公开课、慕课等现代信息技术教学手段的应用, 教师必须改变传统的授课模式, 学会视频制作技巧, 做优教学设计。通过微课学习既满足了学生自主学习的兴趣, 也让教师在制作微课过程中不断提高课堂效率, 还可以让学生在课下依照自己的步调观看课程视频, 通过在线论坛与教师和同学进行交流沟通。在课堂教学环节, 学生在教师的指导下对疑难问题进行分组对话和讨论, 既而问题解决、实践应用和提问与答疑, 这样将传统的课堂进行翻转, 真正做到了因材施教。
3. 加强协作学习。
教学模式的改变, 需要学生积极参与团队的共同任务来促进学生对知识的理解与掌握, 这将提高学生的学习能力、表达能力、思辨能力以及团队协作能力。
混合式学习目前已成为国内外高校教学应用的重要形式。相信在《非线性编辑》课程教学中通过增加交互式视频课件、微课和优化Moodle网络平台等措施, 能真正实现混合式学习, 提高学生的自主学习和协作学习能力, 使教学达到更优化。
参考文献
[1]Graham, C.R.Blended learning systems:definition, current trends, and future directions.In Handbook of Blended Learning:Global Perspectives, Local Designs[M], edited by C.J.Bonk and C.R.Graham, San Francisco, CA:Pfeiffer Publishing, 2006:3-21.
[2]Elizabeth Stacey, Philippa Gerbic.Success factors for blended learning[DB/OL].2010.http://www.ascilite.org.au/conferences/melbourne08/procs/stacey.pdf
[3]百度百科[EB/OL].http://baike.baidu.com/view/283396.htm.
[4]赵国栋, 李志刚.混合式教学与交互式视频课件设计教程[M].北京:高等教育出版社, 2013:4.
[5]赵国栋.信息技术究竟给教学带来什么.2013北京大学“新媒体与学习”暑期学校[Z].北京, 2013.
非线性信号模型概述 篇9
我们生活在信息时代。每个中国人平均每天要接触到35GB以上有价值的信息。这些信息来自于不同的传输媒介, 包括电视、广播、报纸、杂志和万维网, 采取了各种各样的形式, 例如音频、图像、视频、文本等。可靠的有意义的信息的增加是现在社会最大的推动力, 在不同的领域都产生了戏剧性的有形的影响, 如物理、医药、金融、娱乐等。
呈指数增长的海量可用数据带来了一个信息泛滥的问题。一个例子是, 一个普通的石油和天然气公司每天产生10GB以上的数据, 这些数据的来源是他们部署在几千个油井上的数字传感器。再举一个例子, 伦敦的200交通摄像头每天生成并且向中央监控单元发送超过8 TB的视频数据[1]。这种程度的信息量管理给科学和工程带来了重大的挑战。
信号处理可以被非常笼统地称为信息的科学[2]。信号处理研究在设计和分析上的核心部分是:高效记录物理现象的传感器;将数据转化为计算机可以识别的形式的数值方法;允许有价值数据操纵和传输的系统;将处理过的数据转化为对应的有意义的信息的算法。然而, 目前的信息系统为了满足现今的需求所产生的海量数据面临着巨大的压力。传统的信号处理方法并不总是有足够的能力去面对现今的挑战。所以, 信号处理的研究正处于一个十字路口。
2 信号模型
我们现在关注的一个概念叫信号模型。通俗地说, 一个模型是对一类我们感兴趣的给定信号的不同点的数学描述。根据信号大小的不同, 我们感兴趣的信号通常可以通过有限的几个自由度或者有物理意义的参数来表示。换句话说, 我们感兴趣的信号通常展现出简洁而又结构化的一面。
信号模型是一种包含我们感兴趣信号简洁的性质的正式的数学结构。核心前提是, 这样一个简洁的数学描述可以使我们更好地理解问题, 然后给出有效的数值算法结果。信号模型的概念构成了所有信号处理方法的基石, 有时这种现象非常巧妙。例如, 著名的奈奎斯特采样定理[3]指出, 一个有限带宽的连续信号在采样率是带宽两倍的情况下可以被精确地重构。这里最重要的假设就是我们感兴趣的信号是带限的, 因此, 所有的信号模型都是基于带限信号的。
绝大多数传统信号和图像处理, 通信和控制技术都是基于线性模型的。这种模型反映的最基本的概念是, 两个信号线性叠加之后产生的信号同样是在物理上有意义的。线性模型可以解释和分析很多核心的信号处理原理。例如, 支持奈奎斯特采样定理的带限假设, 实际上隐含反映出了一个线性模型。这从带宽为W的信号的线性叠加的结果的带宽依然是W可以很显然地反映出来。
从数学的角度来看, 一个线性信号的几何概念可以很简单地描述为, 一个包含信号的高维信号空间的子空间。传统的信号处理技术, 如信号采样, 重建, 滤波和处理, 可以通过线性子空间的算子进行分析。基于子空间的方法也被成功地应用到其他方面的信号处理应用中, 如估计、跟踪、检测和去噪。线性建模方法不仅使信号能够被分析, 还同时对系统的设计和实施产生了影响。大型的复杂的数字和模拟信号处理系统通常使用简单的模块叠加起来, 而这些模块通常都是线性的。此外, 一些重要的控制系统通常被设计为遵循线性方案。
3 非线性模型
在一些应用领域, 线性假设是个性能非常差的近似假设。一个简单的例子就是现实世界中的图像处理, 两个真实世界图像的线性叠加一般情况下无法构成一幅有意义的图像。在这种情况下, 使用传统的方法比如线性滤波器对噪声图像进行处理经常会产生大量的模糊和其他的人为影响。所以, 对立与传统线性方法的思路, 发展一种非线性的方法和模型并且应用于设计和分析信号处理系统是非常有必要的。
过去很久时间内, 核心信号处理原理的发展是不包括非线性信号模型的。在很多典型的问题中, 通常还是要引入一个线性假设, 即使这个假设的性能比较差。非线性模型和方法通常被用在一些大型公司验证他们的数值计算能力, 因为非线性算法通常需要极其巨大的计算量, 此外, 非线性计算的性能很难给出一个分析上的标准。
但是有点让人吃惊的是, 在21世纪的前十年里, 很多一般性的非线性模型的研究工作取得了很重要的成果, 特别是在机器学习、应用数学等领域。非线性模型带来的灵感让信号处理系统的设计和分析产生了根本的变化。非线性模型应用的核心问题和应用线性模型的传统信号处理原理有着显著的不同。我们这里提到三种不同的比较广为人知的非线性模型的例子。
3.1 稀疏信号模型
稀疏信号模型在一些信号处理应用中比较流行。在这里, 一个基本的假设就是我们感兴趣的信号是稀疏的, 在一个固定的正交基底或者一个冗余字典下可以表示为很少的元素的和。但是, 稀疏信号的线性叠加不一定是稀疏的, 这说明了稀疏信号模型的确是一种非线性模型。实际上, 可以证明稀疏信号存在于整体信号空间的某些高度非线性子空间的并集中。稀疏信号模型已经被应用很多领域比如医学成像, 高速通信, 雷达, 地理信息系统等。压缩感知[4]作为一种新的处理稀疏信号的理论, 给稀疏信号处理带来了极大的发展。
3.2 低秩矩阵模型
低秩矩阵模型在信号处理和机器学习中比较引人注目。在这种模型中, 我们可以使用低秩矩阵和张量来对模型中的数据进行很好的近似。和稀疏信号一样, 低秩矩阵属于所有矩阵的环境向量空间的高度非线性子集。实际上, 低秩矩阵可以看做是一个多项式方程的一个解, 从代数上而言, 就是一个矩阵空间的子空间。低秩矩阵模型主要被用于多目标学习, 视频处理和系统识别领域中。
3.3 流形模型
流形模型已经引起了广泛的关注, 特别是在图像处理领域。在这里, 一个基本假设是接受到的高维的数据是一个仅有少数参数控制的光滑函数。利用这个简洁的性质, 我们甚至有可能打破“高维诅咒”, 而“高维诅咒”是现实世界机器学习中的一个无法解决的问题。然而这个泛函表示和它的相关参数也具有潜在的高度非线性。所以, 从几何的观念看, 数据是被表示在一个高维空间的非线性的低维子流形上的。利用流形结构可以解决大量的诸如机器学习, 计算机视觉, 感官数据融合等问题。
4 总结与展望
基于线性模型的信号处理已经发展到了瓶颈, 但是基于非线性模型的信号处理比如压缩感知近年来得到了长足的发展。在不远的将来, 随着计算机的进步, 非线性信号处理必然会得到更加广泛的应用。
参考文献
[1]A.Szalay and J.Gray, “2020 computing:Science in an exponential world, ”Nature, vol.440, no.7083, pp.413-414, 2006.
[2]C.Shannon, “Communication in the presence of noise, ”Proc.Inst.Radio.Eng., vol.37, no.1, pp.10-21, 1949.
[3]H.Nyquist, “Certain topics in telegraph transmission theory, ”Trans.AIEEE, vol.47, no.2, pp.617-644, 1928.
【非线性混合模型】推荐阅读：
非线性混合效应模型法05-13
混合线性模型07-25
线性混合随机效应模型05-14
非线性函数模型07-29
非线性组合预测模型08-21
非线性有限元模型07-23
混合模型05-12
用MATLAB回归非线性模型参数10-02
混合多模型08-13
有限混合模型09-17