分子量分解问题.doc

分子量分解问题
摘要
本文是针对拥有计算机或不拥有计算机的情况下，求解已知蛋白质分子量的氨基酸组成。为此，我们分别在理论基础上建立了n元一次线性方程模型一和根据实际情况分别补充两个不同的实际约束条件的简化模型二和模型三。
对于模型一即n元一次线性方程模型，在求解过程中，考虑到采用穷举法时计算步骤冗杂，且时间复杂度随X值的增大而呈指数增长，难以较快速度求得数目庞大的解。所以我们另辟蹊径，用C语言编写递推算法，建立可查蛋白质分子量X是否可被第i种氨基酸分子量a[i]拆分的表格，创建查表法，将计算机时间复杂度由指数形式降为多项式形式，提高了运行效率，尤其是针对不拥有计算机的情况，实验人员根据该表格查阅，可提高对蛋白质分子量的拆分速度。相应地，在拥有计算机的情况下，在查表法的基础上继续编写算法，最终实现对X的分解。以X=950为例，共求得 16885个解。
尽管采用查表法可相对提高效率，但模型一的局限在于求得的解数目庞大，不切合实际，实际意义不大。为此我们根据现有的科学数据知识，补充实际情况下的约束条件，在模型一基础上，分别以已知蛋白质中氮含量范围和准确值为约束条件建立简化模型二和模型三，这样大大缩小了解的数量范围，甚至可求得精确解，符合科学实际，使问题得到合理的解决。
关键字 查表法 n元一次线性方程 穷举法 递推算法 简化模型
分子量分解问题
一、 问题重述
生命蛋白质是由若干种氨基酸经不同的方式组合而成。在实验中，为了分析某个生命蛋白质的分子组成，通常用质谱实验测定其分子量x (正整数)，然后将分子量x分解为n个已知分子量a[i](i=1,.......,n)氨基酸的和的形式。某实验室所研究的问题中：
n=18, x 1000
a[i](i=1,.......,18)分别为57, 71, 87, 97, 99, 101, 103, 113, 114, 115, 128, 129, 131, 137, 147, 156, 163, 186
要求针对该实验室拥有或不拥有计算机的情况作出解答。
二、 问题分析
本文重在解决将蛋白质分子如何分解为18种氨基酸的问题。考虑到若使用穷举法对蛋白质分子量进行穷举，计算机运行时间长，且随着分子量x的不断增大，时间复杂度指数增长，尤其是在要求实验室不拥有计算机的情况下，繁杂的数学求解过程使得手工运算量冗杂庞大。为此我们力求寻找其他的方法模型以减少运算工作量，降低问题的复杂度。
在实验室的实际环境下，考虑实验室实际具体情况和能够获得的相关专业实验数据，可进一步获得条件，并结合实验数据在原有基础上构建简化模型求解。
三、 模型假设
1、 所给氨基酸和蛋白质的分子量a[i]和x的值都是精确的；
2、 忽略氨基酸形成蛋白质的过程中脱掉水分子的质量的影响；
3、 各种氨基酸之间没有相互影响，即不存在某种氨基酸的存在已另一种氨基酸的存在为前提的额情况；
4、 在蛋白质中，每种氨基酸出现的概率相等。
四、 符号系统
a[i] 第i种氨基酸的分子量
yi 组成蛋白质的第i种氨基酸的数量
b[i] 第i种氨基酸中含有的氮原子个数
X 已知的某个生命蛋白质的分子量
Pn 已测定的蛋白质的含氮量
五、 模型建立
（一）模型一 n元一次线性方程模型
5.1．1 模型的建立
由蛋白质分子量X和氨基酸分子量a[i]测定蛋白质的组成，假设yi为第i种氨基酸的数量，实际为求解i =X的正整数解yi。故在无实际约束条件的情况细心啊，建立n元一次线性模型：
i =X
yi>=0，yi整数

5.1.2 穷举法求解
在拥有计算机的情况下，通常地，求解该模型采用穷举法，用C语言编写循环算法求解。但是无论采用穷举法循环思想还是动态规划的算法编写程序，求解过程中始终时间复杂度大，耗费时间长，当X为很大的数值时，该方程的解的个数随X的增大以指数的倍数增长，相应的，计算机求解难以在较快速度下得出全部整数解。在尝试过程中，由以下表格数据可以看出：

X（分子量）
Yi（解的个数）
Z（运行的时间/秒）
100
0
0.00033632
150
0
0.00087844
200
4
0.0022
250
3
0.0077
300
14
0.0168
350
7
0.0359
400
45
0.0727
450
44
0.1367
500
158
0.276
550
202
0.4814
600
522
0.9853
650
694
1.7929
700
1508
3.0323
750
2197
4.3742
800
4291
7.0118
850
6337
9.8778
900
11249
17.4341
950
16885
25.8158
1000
28268
42.884

表一：用穷举法求得蛋白质分子量、解的个数和程序运行时间对应表

并且根据查阅资料显示，蛋白质分子量通常在几千甚至几万以上道尔顿，那么利用穷举法对该模型进行求解显然耗时耗力，效率低下，且不切实际，实际意义不大，为此我们放弃该方法的使用。
5.1.3 查表法
1.首先考虑在实验室不拥有计算机的情况下。
从理论角度分析，若人工穷举18种氨基酸的全部组成情况不仅运算量和耗时巨大，而且不切实际。为此，我们考虑创建一种工具表格，通过相关计算和C语言编程，在表格中提供分子量x从1到1000的值，找出当前的x值可分解成18种氨基酸中的哪些氨基酸，其中，1表示X可由a[i]分解，0表示X不可由a[i]分解。实验人员按照这种方法操作查阅，能够以相对较快甚至最快的速度找出所有解的情况。表格模型建立方法如下：
为测定分子量X的蛋白质的氨基酸组成，我们可以将其抽象为对X数值进行拆分，最终将X拆分为几个整数（氨基酸分子量）和的形式。此处，将X可拆分定义为X最多只能由a[i=]7, 71, 87, 97, 99, 101, 103, 113, 114, 115, 128, 129, 131, 137, 147, 156, 163, 186这18个整数组成，也就是说，一个数X如果可以拆分，那么这个数就可以拆分成一个可以继续拆分的数与a中的一个元素的和，继续拆分下去，便可以将一个数拆分成a中的各个元素之和。为避免拆分过程中出现重复拆分，我们对拆分的顺序加以限制，使拆分过程中最先前的元素总不大于后一个元素。
首先根据数据生成两个矩阵B和C，B矩阵是用来存储1到X值的n行1列矩阵，C矩阵是n行18列矩阵，均用来判断1~X这些数是否可以分解。若对于任意的一个数n（1<=n<=x）,n如果可以分解，则B(1,n)=1;否则B(1,n)=0;
与上面思考的过程相反，在C语言编写过程中采用递推的过程：第一步，通过外层设置循环变量p从1到18的循环，使A(p)呈递增趋势；
第二步，令循环变量q从1到x/A(p)循环,使在在(x*B(x)+ q*A(p))<x的情况下，赋值B(x*B(x)+ q*A(p)) = 1， C(x*B(x)+ q*A(p),p) = 1；
第三步，通过18层的循环，将1到1000中的所有可以分解的数进行了一次遍历。结果B中的所有可以分解的数值都被标记为1；
对于C(i,j)中的第j列若为1则表示，i可以分解成a(j)与一个可以拆分的数的和，并且继续拆分下去的时候，其组成元素的值不超过a(j)。（具体代码见附录）
A=[57,71,87,97,99,101,103,113,114,115,128,129,131,137,147,156,163,186];
temp = 1000;
B = zeros(1,temp);
C = zeros(temp,18);
for p = 1:18
for q = 1: floor( temp/A(p) )
for x = 1:temp
if (x*B(x)+ q*A(p)) <= 1000
B(x*B(x)+q*A(p))=1;
C(x*B(x)+q*A(p),p)=1;
end
end
end
end
（所建表格见附件）
具体查表方法如下（以x=200为例）：
首先由表中数据得，a[6],a[7],a[8],a[12]=1，表明该蛋白质可含有第6种，第7种，第8种和第12种氨基酸。在确定含有第6种，分子量为101的氨基酸的情况下，下一步查找x’=200-101=99的分解情况，由表查得x’=99时有且仅有a[5]=1,即该蛋白质还含有分子量为99的第五种氨基酸，故得出该氨基酸的一种分解情况为200=101+99.
在含有分子量为103，第7种氨基酸的情况下，x’=200-103=97,由表查得x’=97时有且仅有a[4]=1,即该蛋白质还含有分子量为97的第四种元素，该蛋白质的第二种复分解情况为200=103+97.同理可求得x=200的其他分解为113+87,129+71.
所以，分子量为200的蛋白质共有4种解。（如图所示）
X=200

a[7]
103
a[8]
113
a[12]
129
a[6]
101

a[2]
71
a[3]
87
a[4]
97
a[5]
99

2、在实验室拥有计算机的情况下
利用先前已建立的表格算法，根据人工查表法利用matlab编写查表算法，可以较快速度得到所求所有解。以X=1000为例，共求得28268个解。
X（分子量）
Yi（解的个数）
Z（运行时间）
100
0
0.021067
150
0
0.019902
200
4
0.025353
250
3
0.023328
300
14
0.023049
350
7
0.025171
400
45
0.028511
450
44
0.028344
500
158
0.030984
550
202
0.035563
600
522
0.054593
650
694
0.066183
700
1508
0.125385
750
2197
0.240396
800
4291
0.524461
850
6337
0.881089
900
11249
2.108896
950
16885
3.994522
1000
28268
9.817489

表二 用查表法求得蛋白质分子量、解的个数和程序运行时间对应表

图一 查表法与穷举法的时间随X值的比较图
由表一和表二，图一和图二可以清晰看出，X值很小时穷举法与查表法运行时间都在1秒之内，然而随着X值的增大，从X=600后，穷举法所用的时间开始陡增，且运行时间远远大于查表法。为此，对于取X较大的值，采用查表法运行时间更短，分解速度更快，效率更高，符合实际的需要。
（二）补充约束条件的简化模型
第一种通用模型由于无实际约束条件所以解的数量巨大，人工求解工作量耗时费力，实际意义不大。为此我们考虑是否可以根据实验室的实际情况补充相应的约束条件，以此剔除不满足实际情况的解，减少求解数目，简化求解工作量，加快求解速度。
经过查阅相关资料，根据现有实验室的实际具体条件设备，可获知蛋白质的氮含量。
5.2.1模型二 已知蛋白质的氮含量范围
根据现有的科学知识知道“生命蛋白质中氮含量约在15%~17%之间”，假设b[i]为第i种氨基酸中含有氮原子的数量，据此建立含氮量模型：
i=X
0.15<=<=0.17
在基于原模型查表法的基础上，根据补充的含氮量的约束条件，用C语言编程，在查表法得出的所有解中逐一剔除不满足含氮量的解，最终利用matlab得出所有解得情况如下：

X
Yi解的个数
100
0
150
0
200
0
250
1
300
0
350
3
400
0
450
9
500
115
550
104
600
56
650
95
700
732
750
341
800
589
850
2317
900
3839
950
2470
1000
13421

5.2.2模型三 已知蛋白质中氮元素的准确含量
经过查阅资料，我们了解到在生物实验室中可以通过凯氏定氮法、显色反应等科学实验比较精确地测定蛋白质的含氮量Pn.
i=X
=Pn
当我们获得蛋白质氮元素的精确含量后，通过该简化模型运用查表法可大大减小解的范围，并且得到精确解。
六、 模型分析
本文提出了一种无实际约束条件的n元一次线性模型，并在原有模型基础上根据科学数据补充实际约束条件，进一步提出简化模型
对于模型一，在纯理论的条件下提出，缺少实际约束条件，求得的解数量庞大，无论对于拥有计算机还是没有计算机的情况，求解过程中时间复杂度高，并且随着X值的不断增大，求得的解得数目庞大，导致求解效率低下，不符合实验室实际具体情况，实际意义不大。
对于模型二，根据实验室现有的实验测试仪器设备和相关科学知识数据，我们在补充实际约束解后，缩小了解的范围，使得求解更为精确，将理论上升为实际的需要，满足实际的要求。
为此，选择模型二更贴合实际，满足实际的需要。
在求解模型的方法上，利用穷举法计算过程冗杂，时间复杂度随X的值增大呈指数增长，尤其对于不拥有计算机的情况，人工求解工作量庞大，耗时巨大，工作效率低下。
而查表法根据之前制定的数据表格，能够有效地查取可以拆分的X值，实现对X的分解，尤其对于没有计算机的情况，对X的分解次数对应X的所有解得情况，可最大程度上减小工作量，加快求解速度，提高效率。
七、 模型推广
此类模型还可以运用在投资、运输、订购物资等方面运用。在投资方面，假设拥有一定资金可以运转，可以购买不同公司的股票，通过本模型可以计算出所有可以投资的种类，然后通过对每种股票风险的计算得出最佳投资方式；在运输方面，有不同种运输方式货运输线路，在资金一定的情况下计算出以何种组合的方式或路线最合适，然后再反过来计算最小资金；在订购物资方面，有不同厂家提供不同价格不同质量的物资（但都符合订购的最基本要求），可以通过本模型的计算得出最合理有又最省资金的订购方式。
八、 结论
本文是针对拥有计算机或不拥有计算机的情况下，求解已知蛋白质分子量的氨基酸组成。为此，我们分别在理论基础上建立了n元一次线性方程模型一和根据实际情况分别补充两个不同的实际约束条件的简化模型二和模型三。
在模型的求解过程中，由于穷举法的计算过程冗杂，且程序时间复杂度随着X值的增大而呈指数增长，求解运行时间长，效率低下。为避开这些弊端，我们另辟蹊径地采用查表法，建立可查蛋白质分子量X是否可被第i种氨基酸分子量a[i]拆分的表格，相比于传统的穷举法，我们将计算机时间复杂度由指数形式降为多项式形式，时间复杂度不会随着X增大而陡增，提高了效率， 这是查表法的优点所在。运用模型一求得结果，以X的值50为梯度，解的个数为：
X（分子量）
Yi（解的个数）
100
0
150
0
200
4
250
3
300
14
350
7
400
45
450
44
500
158
550
202
600
522
650
694
700
1508
750
2197
800
4291
850
6337
900
11249
950
16885
1000
28268
但模型一的局限在于仅在理论的前提下建立，因此求得的可行解数目庞大，而且不切合实际，实际意义不大。为此，补充约束条件已知蛋白质的氮含量范围的模型二和以已知蛋白质中氮元素的准确含量为约束条件的模型三的建立大大缩小了解的数量范围，甚至可求得精确解，符合科学实际，使问题得到合理的解决。
以下是模型二求解的情况：
X
Yi解的个数
100
0
150
0
200
0
250
1
300
0
350
3
400
0
450
9
500
115
550
104
600
56
650
95
700
732
750
341
800
589
850
2317
900
3839
950
2470
1000
13421

九、 参考文献
[1] 程龙，张云军，赵蕊 蛋白质氨基酸的组合问题 百度文库
[2] 王晓东 计算机算法设计与分析 电子工业出版社

附录
附录一 生成表格matlab源代码
tic%计算时间和toc对应
clc%清屏
clear%删除数据
A=[57,71,87,97,99,101,103,113,114,115,128,129,131,137,147,156,163,186];
temp = 1000;
B = zeros(1,temp);
C = zeros(temp,18);
for p = 1:18
for q = 1: floor( temp/A(p) )
for x = 1:temp
if (x*B(x)+ q*A(p)) <= 1000
B(x*B(x)+q*A(p))=1;
C(x*B(x)+q*A(p),p)=1;
end
end
end
end
in = 777;
flag = 0;
if B(in) == 0
disp('无法分解');
else
disp('可以分解');
flag = 1;
end
toc
附录二 穷举法源代码
#include<time.h>
int a[18]={0};
int b[18]={0};
int c[]={57,71,87,97,99,101,103,113,114,115,128,129,131,137,147,156,163,186};
int sum=0;
FILE *fp2;
void search(int m,int k)
{
int i;
int t;
int temp;
if(m==0)
{
for(i=0;i<18;i++)
{
fprintf(fp2,"%5d",a[i]);
}
fprintf(fp2,"\n");
sum++;
}
else if( (m > 0) &&(k<18))
{
t=0;
temp = m/c[k];
for(a[k] = 0;a[k]<=temp;a[k]++)
{
search(m-a[k]*c[k],k+1);
}

}
a[k]=0;
}
void main()
{
fp2 = fopen("result.txt","w");
int n=1000;
double clock_t,start,finish;
double totaltime;
start=clock();
search(n,0);
fprintf(fp2,"%d\n",sum);
finish=clock();
totaltime=(double)(finish-start)/CLOCKS_PER_SEC;
printf("%f\n",totaltime);

}
附录三 拥有计算机情况下，模型一中查表法分解蛋白质分子量X源代码：
%分子量分解
tic%计算程序运行时间
clc
clear
a=[57,71,87,97,99,101,103,113,114,115,128,129,131,137,147,156,163,186];
temp=1000;%分子量上限
B=zeros(1,temp);%若B(x)=1，则代表分子量为x时有解，初始化为全0
C=zeros(temp,18);%记录最后所加分子量的值，若C(x，y)=1，则代表分子量为x时，其中一种分解方式最后所加分子量为y，初始化为全0
for k=1:18 %B，C矩阵运算开始
for j=1:fix(temp/a(k))
for i=1:temp

if i*B(i)+j*a(k)<=temp
B(i*B(i)+j*a(k))=1;
C(i*B(i)+j*a(k),k)=1;

end
end
end
end
%B，C矩阵运算结束
sr=1000;%假设所测分子量为999
if ~B(sr) %B(sr)=0，代表无法分解
disp('无解')
zgs=0
else %可以分解
sr1=sr; %初始化sr1
temp2=zeros(1,19);
temp2(1,19)=18;
zgs=1;%分解方案个数初始化为1
for j=1:17 %分解步骤最多为1000/57=17.54,故大循环次数取17
%计算分解方案开始
zszgs=0;
for k=1:zgs

sr1=sr-temp2(k,1:18)*a';

if sr1
gs=sum(C(sr1,1:temp2(k,19)));

temp1=ones(gs,19);
for l=1:gs
temp1(l,:)=temp2(k,:);
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

js=0;
for i=1:temp2(k,19)
if C(sr1,i)
js=js+1;
temp1(js,i)=temp1(js,i)+1;
temp1(js,19)=i;
end
end
temp3(zszgs+1:zszgs+gs,:)=temp1;
else
gs=1;
temp3(zszgs+1,:)=temp2(k,:);

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
zszgs=zszgs+gs;
end

temp2=temp3;
zgs=zszgs;
end
end
%计算分解方案结束

disp('输入的分子量')
sr

disp('分解结果')
%temp2(:,1:18)
% disp('验证')
% temp2(:,1:18)*a'
disp('解的个数')
zgs
disp('运算时间')
toc
附录四 模型二查表法分解蛋白质源代码
tic
clc
clear
A=[57,71,87,97,99,101,103,113,114,115,128,129,131,137,147,156,163,186];
temp = 1000;
B = zeros(1,temp);
C = zeros(temp,18);
D = [1,1,1,1,1,1,1,1,2,1,2,1,1,3,2,2,1,1];
for p = 1:18
for q = 1: floor( temp/A(p) )
for x = 1:temp
if (x*B(x)+ q*A(p)) <= 1000
B(x*B(x)+q*A(p))=1;
C(x*B(x)+q*A(p),p)=1;
end
end
end
end
sr=200;%假设所测分子量为999
if ~B(sr) %B(sr)=0，代表无法分解
disp('无解')
else %可以分解
sr1=sr; %初始化sr1
temp2=zeros(1,19);
temp2(1,19)=18;
zgs=1;%分解方案个数初始化为1
for j=1:17 %分解步骤最多为1000/57=17.54,故大循环次数取17
%计算分解方案开始
zszgs=0;
for k=1:zgs

sr1=sr-temp2(k,1:18)*A';

if sr1
gs=sum(C(sr1,1:temp2(k,19)));

temp1=ones(gs,19);
for l=1:gs
temp1(l,:)=temp2(k,:);
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

js=0;
for i=1:temp2(k,19)
if C(sr1,i)
js=js+1;
temp1(js,i)=temp1(js,i)+1;
temp1(js,19)=i;
end
end
temp3(zszgs+1:zszgs+gs,:)=temp1;
else
gs=1;
temp3(zszgs+1,:)=temp2(k,:);

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
zszgs=zszgs+gs;
end

temp2=temp3;
zgs=zszgs;
end
end
temp4 = zeros(1,19);
%计算分解方案结束
k = 1;
for p = 1:zgs
tt=temp2(p,1:18)*D'*14;
if (tt/sr)>0.15 &&(tt/sr)<0.17
temp4(k,:)=temp2(p,:);
k = k + 1;
end
end

disp('输入的分子量')
sr

disp('分解结果')
% xlswrite('表格.xls',C(:,1:18))
% xlswrite('result.xls',temp4(:,1:18))
% temp4(:,1:18)
% disp('验证')
% temp2(:,1:18)*a'
disp('解的个数')
k
disp('运算时间')
toc