刘连香 摘 要：本文介绍了数学建模和数据曲线拟合的原理方法，并给出了非线性曲线拟合的最小二乘法的MATLAB实现过程，解决了一个医学上的实际问题――给药方案制定问题。
关键词：曲线拟合； 最小二乘；给药方案
Abstract:This paper introduces the mathematical model and data fitting principle method, and gives the nonlinear curve fitting by least square method in the process of MATLAB implementation, solves a practical medical problems - delivery plan problem.
Key words:curve fitting; least squares; dosing regimens
在很多科学实验、工程数学中，常常需要确定一个变量依存于另一个或更多的变量的关系，即函数。但实际上确定函数的形式（线性形式、乘法形式、幂指形式或其它形式）时往往没有先验的依据。只能在收集的实际数据的基础上对若干合乎理论的形式进行试验，从中选择一个最能拟合有关数据（即最有可能反映实际问题）的函数形式。下面介绍用MATLAB求解曲线拟合的实际问题。
一、曲线拟合的最小二乘法
基本思路：对科学实验有一组实验数据(xi,yi),i=1,2,…,m ，要求一个函数f(x)=a1f1(x)+a2f2(x)+…+amfm(x)与所给数据(xi,yi),i=1,2,…,m 拟合，使函数与所有数据点最为接近，其中fk(x)是事先选定的一组函数，ak是待定系数（k=1,2…,m,m 用最小二乘法求拟合曲线时，首先要确定f(x)的形式，这不单纯是数学问题，还与所研究问题的运动规律及所得观测数据有关。通常要从问题的运动规律及给定数据描图，确定f(x)的形式，并通过实际计算选出较好的结果。
二、MATLAB求解过程
在最小二乘法拟合中，若要寻求的函数f(x)是任意的非线性函数，则称为非线性最小二乘拟合。MATLAB6?3的优化工具箱中提供了求非线性最小二乘拟合函数：lsqcurvefit。使用这个命令时，要先建立M-文件fun?m，在其中定义函数F(x)。
设已经xdata=(xdata1,xdata2,…,xdatan)，ydata=(ydata1,ydata2,…,ydatan)，lsqcurvefit用以求参量x（向量）的向量值函数：
F(x,xdata)=(F(x,xdata1),…,F(x,xdatan))T中的参变量x（向量），使得
∑ni=1(F(x,xdatai)－ydata)2最小
输入格式为：x=lsqcurvefit(‘fun’,x0,xdata,ydata)
输出目标函数值格式：f=fun(x,xdata)，其中x0为迭代的初值。
实 例 给药方案制定问题：一种新药用于临床之前，必须设计给药方案。在快速静脉注射的给药方式下，所谓给药方案是指：每次注射计量多大，间隔时间多长。药物进入机体后血液输送到全身，在这个全过程中不断地被吸收、分布、代谢，最终排出体外，药物在血液中的浓度，即单位体积血液中的药物含量，称为血药浓度。在最简单的一室模型中，将整个机体看作一个房室，称中心室，室内血药浓度是均匀的。快速静脉注射后，浓度立即上升；然后迅速下降；当浓度太低时，达不到预期的治疗效果；当浓度太高，又可能导致药物中毒或副作用太大。临床上，每种药物有一个最小有效浓度c1和一个最大有效浓度c2.设计给药方案时，要使血液浓度保持在c1-c2之间。本题c1=10，c2=25(ug/ml)。
通过实验，对某人用快速静脉注射方式一次注入该药物300mg后，在一定时刻t(h)采集血药，测得血药浓度c（ug/ml）如下表：
解： 问题分析：要设计给药方案，必须知道给药后血药浓度随时间变化的规律。
在半对数坐标中做出t与c的关系图（如右下图），此图说明血药浓度数据符合负指数变化规律。
1.求血药浓度变化规律
假设：（1）机体看作一个房室，室内血药浓度均匀――一室模型。
（2）药物排除速率与血药浓度成正比，比例系数k（>0）。
（3）血液容积为v，t=0时注射剂量为d，则t=0血药浓度为d/v。
由假设（2）得： dcdt=－kc
由假设（3）得：c(0)=d/v
解微分方程：dcdt=－kc
c(0)=d/v
得：c(t)=dve－kt
在此，d=300mg,t及c(t)在某些点处的值见前表，下面通过拟合求出参数k，v用非线性最小二乘拟合c(t)--用lsqcurvefit函数。
（1）先用M-文件curvefun3?m定义函数
Function f=curvefun3(x,tdata)
d=300
f=(x(1)\\d)*exp(-x(2)*tdata) %x(1)=v;x(2)=k
(2)再编写主程序lihe2?m
Clear
tdata=［0?25 0?5 1 1?5 2 3 4 6 8］；
cdata=［19?21 18?15 15?36 14?10 12?89 9?32 7?45 5?24 3?01］;
x0=［10,0?5］;
x=lsqcurvefit(‘curvefun3’,x0,tdata,cdata);
f=curvefun3(x,tdata)
x
计算结果：k=0?2420(1/h),v=14?8212(L).
2.制定给药方案
假设：(1)设每次注射剂量D，间隔时间。(2)血药浓度c(t)应c1≤c(t)≤c2。(3)初次剂量D0应加大。给药方案记为：D0,D,τ,则:
(1) D0=vc2,D=v(c2－c1)
(2) c1=c2e－kτ?τ=1klnc2c1
将c1=10, c2=25, k=0?2347, v=15?02代人（1）（2）计算得：
D0=375?5,D=225?3,τ=3?9
故可制定给药方案：D0=375(mg),D=225(mg),τ=4(h)即首次注射375mg，其余每次注射225mg，注射的间隔时间 为4h.（作者单位：江西城市职业学院）
参考文献
［1］ 张韵华,奚梅成,陈效群.数值计算方法和算法［M］.北京:科学出版社,2002。
［2］ 于润伟. MATLAB基础及应用［M］.北京：机械工业出版社, 2003。
［3］ 李庆扬 王能超 易大义 数值分析（第4版） 清华大学出版社，2005。
［4］ 萧树铁 ，数学实验，高等教育出版社，1999。