数学建模新止痛剂生效时间预测模型.doc

数学建模论文

题目：新止痛剂生效时间预测模型

班级：计算机科学与技术1401班

学号：201409824 姓名：马元凯

学号：201409822 姓名：李炳毅

日期：2016年12月15日

目录

摘要 1
一、 问题描述 1
二、 问题分析 2
三、 模型假设 3
四、 模型建立 3
五、 模型求解 6
六、 模型分析 6
七、 模型改进 9
八、 模型评注与推广 11
九、 参考文献 12
十、 附录 12

摘要
某医药公司的新药研究部门为了掌握一种新止痛剂的疗效，进行了相应的研究，根据医药公司给出的相关实验数据，通过MATLAB软件进行数据处理工作，根据拟合曲线和散点图建立基本模型（模型（1））.用MATLAB软件解出该基本模型，再利用残差图两次剔除数据并且回归，得到最佳模型（模型（2））,该模型的拟合度为85.14%.拟合度不是很高，针对不同的性别，引入独立变量的交互作用，对模型进行改进，得到两个拟合度较高的模型（4）和（5），针对男性建立的模型，拟合度为90.87%，针对女性建立的模型，拟合度为97.72%.
本文给出的模型较为科学，若实验数据真实可信，则可以进行药物的推广.

关键字:MATLAB软件，残差图，拟合度,统计回归模型

一、 问题描述
一个医药公司的新药研究部门为了掌握一种新止痛剂的疗效，设计了一个药物试验，给患有同种疾病的病人使用这种新止痛剂的以下4个剂量中的某一个：2 g，5 g，7 g和10 g，并记录每个病人病痛明显减轻的时间（以分钟计）. 为了解新药的疗效与病人性别和血压有什么关系，试验过程中研究人员把病人按性别及血压的低、中、高三档平均分配来进行测试. 通过比较每个病人血压的历史数据，从低到高分成3组，分别记作0.25，0.50和0.75. 实验结束后，公司的记录结果见下表1（性别以0表示女，1表示男）.
请你为该公司建立一个数学模型，根据病人用药的剂量、性别和血压组别，预测出服药后病痛明显减轻的时间.

病人序号
病痛减轻时间／min
用药剂量／g
性别
血压组别
1
35
2
0
0.25
2
43
2
0
0.50
3
55
2
0
0.75
4
47
2
1
0.25
5
43
2
1
0.50
6
57
2
1
0.75
7
26
5
0
0.25
8
27
5
0
0.50
9
28
5
0
0.75
10
29
5
1
0.25
11
22
5
1
0.50
12
29
5
1
0.75
13
19
7
0
0.25
14
11
7
0
0.50
15
14
7
0
0.75
16
23
7
1
0.25
17
20
7
1
0.50
18
22
7
1
0.75
19
13
10
0
0.25
20
8
10
0
0.50
21
3
10
0
0.75
22
27
10
1
0.25
23
26
10
1
0.50
24
5
10
1
0.75
表1

二、 问题分析
一般来说，药物的疗效可以直观的用服药后病痛明显减轻的时间来衡量.在新药推广中，医药公司的新药研究部门设计了一种药物给患有同种疾病的病人使用后，根据病人的用药剂量、性别和血压组别，预测病痛减轻时间的多少来预测止痛药的疗效，这是一个统计回归问题，针对这个问题，我们需要给出合理的假设，尤其是变量的选取，进而预测出药物的疗效.
先将因变量分别与变量进行单独分析，得出两者间的大致函数关系，进一步整合这些关系，得出一个因变量与各个变量间的关系函数模型并进行求解，得出最终结论。

三、 模型假设
1、 假设病人只服用了新型止痛药，未服用其它药物.
2、 假设题中给出的实验数据真实可信，误差很小.
3、 假设24名病人都是在服用新型止痛药的人群中随机选取的.
4、 假设病人在实验阶段吃的食物对新型止痛药无影响.
5、 假设模型中出现的符号含义如下表2所示.

符号
含义
单位

用药剂量

性别
女-0，男-1

血压组别
低-0.25，中-0.5，
高-0.75

病痛减轻时间
min
p
概率值

随机误差

回归系数

置信水平

表2

四、 模型建立
为了大致地分析y与，，之间的关系，首先利用表1-1的数据分别作出y对，和的散点图（见图4-1，图4-2和图4-3的圆点）.
如图1为y对的散点图，图2为y对的散点图，图3为y对的散点图.

图1

图2

图3
由上图可知：y对可用二次函数拟合，拟合后如图4所示.

图4
根据对散点图图1和图4的分析可得出y对的二次函数模型 ；根据图2可得出y对的线型模型；根据图3可得出y对的线型模型 ；结合上述的3个模型可建立如下（1）的多元线性回归模型：
（1）
上式右端的，，称为回归变量（自变量），给定，，时，病痛减轻时间y的平均值为；由表1的数据估计，影响y的其他因素作用都包含在随机误差中，如果模型选的合适，应大致服从均值为0的正态分布.

五、 模型求解
直接利用MATLAB统计工具箱中的regress求解，进一步求出回归系数估计值及其置信区间（置信水平=0.05）、检验统计量R2，F，p，s2，详细数据见表3所示.
参数
参数估计值
置信区间

63.1291
[48.7173，77.5409]

-10.2706
[-14.9243，-5.6169]

5.6667
[-0.0213，11.3546]

-1.5000
[-15.4325，12.4325]

0.5111
[0.1319，0.8903]
R2=0.8275 F=22.7903 p=0.0000 s2=44.3109
表3
由表3得出的模型为:

六、 模型分析
表3显示，R2=0.8275 指因变量y（病痛减轻时间）的82.75%可由模型确定，F值远远超过F检验的临界值，p远小于，因而模型（1）从整体来看是可用的. 表3的回归系数给出了模型（1）中参数的估计值，，的置信区间包含零点，因而对这两个系数的解释是不可靠的，所以需要残差分析.
首次回归所得残差图如下图5所示，可以看出，第3个和第24个数据存在异常，剔除，进行第二次回归.

图5
第二次回归所得残差图如下图6所示，可以看出，第5个数据存在异常，剔除，进行第三次回归.

图6

第三次回归所得残差图如下图7所示.

图7
第三次回归的结果如表4所示.

参数
参数估计值
参数置信区间

55.8121
[44.3918 67.2324]

-8.0962
[-12.0605 -4.1320]

7.1311
[2.5539 11.7083]

-6.3868
[-18.5348 5.7612]

0.4025
[0.0869 0.7181]
R2= 0.8514 F= 22.9153 p= 0.0000 s2= 24.1056
表4
由图7可知数据没有异常项，因此模型基本可用，此时得出的最佳模型应为： （2）

易知因变量y（病痛减轻时间）的85.14%可由模型确定.

七、 模型改进
从以上的分析可以看出，模型的拟合度最大为85.14%，拟合度不是很高，需要进行模型的改进；因为服药期间要考虑生理反应，而性别对生理反应有直观的影响，故改进的模型需要对男女分开进行讨论；模型（2）中回归变量和对因变量y的影响是相互独立的.根据直觉和经验可以猜想和之间的相互作用会对y有影响，不妨简单地用和的乘积代表它们的交互作用，于是在模型（1）中增加一项，得到
（3）

1、 针对男性，利用表1的数据估计模型（3）的系数，利用MATLAB得到如下表5的结果.

参数
参数估计值
参数置信区间

49.8088
[24.4805 75.1372]

-7.8431
[-14.4259 -1.2604]

39.0294
[-1.0850 79.1438]

-7.5882
[-13.6016 -1.5748]

0.6667
[0.1895 1.1438]
R2= 0.9087 F= 17.4206 p= 0.0010 s2= 27.4856
表5
根据表5，可得出此时的最佳模型为：
（4）

对表5的数据进行残差分析，得到如下图8的结果.

图8
从图8可以看出数据无异常项，则由此得出的模型（4）基本可信，病痛减轻时间的90.87%可由模型确定.
2、 针对女性，利用表1的数据估计模型（3）的系数，利用MATLAB得到如下表6的结果.

参数
参数估计值
参数置信区间

36.9395
[22.9221,50.9570]

-5.1686
[-8.8117,-1.5255]

48.3235
[26.1230,70.5240]

-7.4706
[-10.7986,-4.1426]

0.3556
[0.0915,0.6196]
R2=0.9772 F=74.8974 p=0.0000 s2= 8.4184
表6
对表6的数据进行残差分析，得到如图9所示的结果,

图9
从图9可以看出，第8项数据存在异常，但影响不大，所以不用剔除数据，模型基本可信.根据表5得出此时的最佳模型为：
（5）
病痛减轻时间的97.72%可由模型确定.通过以上分析可以确定，经过对模型的改进，针对不同的性别，得出的两种模型（4）和（5）的拟合程度都达到了百分之九十以上，达到了模型改进的目的.

八、 模型评注与推广
评注 从这个实例我们可以看出，建立回归模型要根据已知的数据，从常识和经验进行分析，辅以作图，决定最终的函数模型.用MATLAB软件求解后，作统计分析R2的值是对模型的直观评价，决定模型的拟合程度.本次建立的模型，优点在于，由简到繁，先单独考虑变量对结果的影响，再综合考虑，最后引入交互项对模型进行改进，为病痛减轻时间的预测提供了简洁方便的工具，用简单的形式表示出了事物之间复杂的关系；缺点在于，该模型不能实际的去观测数据，而是使用实验给出的数据，存在一定的误差，还有改进模型中针对不同性别得出的模型拟合程度差别较大，不是很理想；其实在模型建立时，还可以引进其他的数据，如心率等因素，使模型更加真实可信.
推广 我们建立的这个预测模型还可以用于软件开发人员的薪金预测或其它的薪资预测模型.

九、 参考文献
[1] 姜启源, 谢金星, 叶俊. 数学模型.第4版[M]. 高等教育出版社, 2011.
[2] 宋来忠，王志明,数学建模与实验,北京：科学出版社,2005
[3] 姜启源, 谢金星, 叶俊. 数学模型 (第四版) 习题参考解答[M]. 高等教育出版社, 2011.（第十章习题参考答案）

十、 附录
1、图1，图2，图3，图4的MATLAB程序：
y对的散点图
>> x1=[2 2 2 2 2 2 5 5 5 5 5 5 7 7 7 7 7 7 10 10 10 10 10 10];
>> y=[35 43 55 47 43 57 26 27 28 29 22 29 19 11 14 23 20 22 13 8 3 27 26 5];
>> scatter(x1,y,'r');
y对的散点图
>> x2=[0 0 0 1 1 1 0 0 0 1 1 1 0 0 0 1 1 1 0 0 0 1 1 1];
>> y=[35 43 55 47 43 57 26 27 28 29 22 29 19 11 14 23 20 22 13 8 3 27 26 5];
>> scatter(x2,y,'r');
y对的散点图
>> x3=[0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75];
>> y=[35 43 55 47 43 57 26 27 28 29 22 29 19 11 14 23 20 22 13 8 3 27 26 5];
>> scatter(x3,y,'r');
y对的拟合曲线
>> x1=[2 2 2 2 2 2 5 5 5 5 5 5 7 7 7 7 7 7 10 10 10 10 10 10];
>> y=[35 43 55 47 43 57 26 27 28 29 22 29 19 11 14 23 20 22 13 8 3 27 26 5];
>> p=polyfit(x1,y,2);
>> x1x1=linspace(min(x1),max(x1));
>> yy=polyval(p,x1x1);
>> plot(x1,y,'o',x1x1,yy);
2、表3的MATLAB程序：
>> x1=[2 2 2 2 2 2 5 5 5 5 5 5 7 7 7 7 7 7 10 10 10 10 10 10];
>> x2=[0 0 0 1 1 1 0 0 0 1 1 1 0 0 0 1 1 1 0 0 0 1 1 1];
>> x3=[0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75];
>> y=[35 43 55 47 43 57 26 27 28 29 22 29 19 11 14 23 20 22 13 8 3 27 26 5];
>> x=[ones(24,1),x1',x2',x3',(x1.^2)'];
>> [b,bint,r,rint,stats]=regress(y',x)
3、 图5，图6，图7的MATLAB程序：
首次回归分析
>> x1=[2 2 2 2 2 2 5 5 5 5 5 5 7 7 7 7 7 7 10 10 10 10 10 10];
>> x2=[0 0 0 1 1 1 0 0 0 1 1 1 0 0 0 1 1 1 0 0 0 1 1 1];
>> x3=[0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75];
>> y=[35 43 55 47 43 57 26 27 28 29 22 29 19 11 14 23 20 22 13 8 3 27 26 5];
>> x=[ones(24,1),x1',x2',x3',(x1.^2)'];
>> [b,bint,r,rint,stats]=regress(y',x,0.05)

>> rcoplot(r,rint)
第二次回归分析
>> x1=[2 2 2 2 2 5 5 5 5 5 5 7 7 7 7 7 7 10 10 10 10 10];
>> x2=[0 0 1 1 1 0 0 0 1 1 1 0 0 0 1 1 1 0 0 0 1 1];
>> x3=[0.25 0.50 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50];
>> y=[35 43 47 43 57 26 27 28 29 22 29 19 11 14 23 20 22 13 8 3 27 26];
>> x=[ones(22,1),x1',x2',x3',(x1.^2)'];
>> [b,bint,r,rint,stats]=regress(y',x,0.05)

>> rcoplot(r,rint)
第三次回归
>> x1=[2 2 2 2 5 5 5 5 5 5 7 7 7 7 7 7 10 10 10 10 10];
>> x2=[0 0 1 1 0 0 0 1 1 1 0 0 0 1 1 1 0 0 0 1 1];
>> x3=[0.25 0.50 0.25 0.50 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50];
>> y=[35 43 47 43 26 27 28 29 22 29 19 11 14 23 20 22 13 8 3 27 26];
>> x=[ones(21,1),x1',x2',x3',(x1.^2)'];
>> [b,bint,r,rint,stats]=regress(y',x,0.05)

>> rcoplot(r,rint)
4、 模型改进后，相关的MATLAB程序
表5，图8的MATLAB程序
>> x1=[2 2 2 5 5 5 7 7 7 10 10 10 ];
>> x3=[0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75];
>> y=[47 43 57 29 22 29 23 20 22 27 26 5];
>> x=[ones(12,1),x1',x3',(x1.*x3)',(x1.^2)'];
>> [b,bint,r,rint,stats]=regress(y',x)

>> rcoplot(r,rint)
表6，图9的MATLAB程序
>> x1=[2 2 2 5 5 5 7 7 7 10 10 10];
>> x3=[0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 0.25 0.50 0.75 ];
>> y=[35 43 55 26 27 28 19 11 14 13 8 3];
>> x=[ones(12,1),x1',x3',(x1.*x3)',(x1.^2)'];
>> [b,bint,r,rint,stats]=regress(y',x)

>> rcoplot(r,rint)

目 录

第一章 总 论 1
一、项目提要 1
二、可行性研究报告编制依据 2
三、综合评价和论证结论 3
四、存在问题与建议 4
第二章 项目背景及必要性 5
一、项目建设背景 5
二、项目区农业产业化经营发展现状 11
三、项目建设的必要性及目的意义 12
第三章 建设条件 15
一、项目区概况 15
二、项目实施的有利条件 17
第四章 建设单位基本情况 19
一、建设单位概况 19
二、研发能力 20
三、财务状况 20
第五章 市场分析与销售方案 21
一、市场分析 21
二、产品生产及销售方案 22
三、销售策略及营销模式 22
四、销售队伍和销售网络建设 23
第六章 项目建设方案 24
一、建设任务和规模 24
二、项目规划和布局 24
三、生产技术方案与工艺流程 25
四、项目建设标准和具体建设内容 26
五、项目实施进度安排 27
第七章 投资估算和资金筹措 28
一、投资估算依据 28
二、项目建设投资估算 28
三、资金来源 29
四、年度投资与资金偿还计划 29
第八章 财务评价 30
一、财务评价的原则 30
二、主要参数的选择 30
三、财务估算 31
四、盈利能力分析 32
五、不确定性分析 33
六、财务评价结论 34
第九章 环境影响评价 35
一、环境影响 35
二、环境保护与治理措施 35
三、环保部门意见 36
第十章 农业产业化经营与农民增收效果评价 37
一、产业化经营 37
二、农民增收 38
三、其它社会影响 38
第十一章 项目组织与管理 40
一、组织机构与职能划分 40
二、项目经营管理模式 42
三、技术培训 42
四、劳动保护与安全卫生 43
第十二章 可行性研究结论与建议 46
一、可行性研究结论 46
二、建议 47