問題一難度較低，主要是找出對生物活性最具有顯著影響的分子描述符，具有影響意味著該分子描述符與生物活性的相關性越高，因此可以構建相關性分析模型來分別計算各分子描述符與生物活性的關系，找出相關性最大的20個變量即可。目前有同學說用皮爾遜系數做的結果并不是很好，這是因為皮爾遜系數主要適用于單調線性變化的數據，而題目中的數據并沒有明顯的單調性質，因此得到的結果可能并不好，這時候建議大家可以試一下灰色關聯分析模型，他是根據因素之間發展趨勢的相似或相異程度，亦即“灰色關聯度”，作為衡量因素間關聯程度的一種方法。

灰色關聯分析代碼如下

％數據替換為自身研究數據，代碼會自動去除量綱clear;clc;load data.mat;r = size(data,1);c = size(data,2);%第一步，對變量進行預處理，消除量綱的影響avg = repmat(mean(data),r,1);data = data./avg;%定義母序列和子序列Y = data(:,1); %母序列X = data(:,2:c); %子序列Y2 = repmat(Y,1,c-1); %把母序列向右復制到c-1列absXi_Y = abs(X-Y2)a = min(min(absXi_Y)) %全局最小值b = max(max(absXi_Y)) %全局最大值ro = 0.5; %分辨系數取0.5gamma = (a+ro*b)./(absXi_Y+ro*b) %計算子序列中各個指標與母序列的關聯系數disp("子序列中各個指標的灰色關聯度分別為：");ans = mean(gamma)

針對問題二：該問題主要是構建生物活性與20個分子變量之間的邏輯關系，主要有兩種思路，第一種是構建多元線性模型，及構建y=a1x1+a2x2+...+a20x20關系，然后求解a1,a2,...a20系數即可；第二種思路是構建神經網絡模型，建立變量與生物活性之間的非線性擬合關系，這里最好別用BP神經網絡，因為輸入數據的參數較多，很容易陷入局部最優解，可采用遺傳算法優化的BP神經網絡模型或粒子群算法優化的神經網絡模型等，建立上述邏輯關系后將test表中數據帶入即可。

有同學說這個關系如何構建，能不能用回歸分析的思路，回歸分析肯定是用線性回歸，非線性對于這么多的指標難度太大而且也不好把握，但如果線性回歸其關系式又太過簡單（畢竟醫學領域用線性關系構建的可能性很小）因此回歸算法在這里并不推薦。

還是主推神經網絡之類的非線性擬合模型，對于未知關系的擬合不借助任何參考關系式，而是采用權值鏈接在一起，精度會更高一些，但注意過擬合和誤差分析問題，并不局限于用BP神經網絡，徑向基網絡，LM－BP神經網絡等也可以試一下

?

%程序一：GA訓練BP權值的主函數
function net=GABPNET(XX,YY)

%--------------------------------------------------------------------------

% ?GABPNET.m

% ?使用遺傳算法對BP網絡權值閾值進行優化，再用BP算法訓練網絡

%--------------------------------------------------------------------------

%數據歸一化預處理

nntwarn off

XX=[1:19;2:20;3:21;4:22]';

YY=[1:4];

XX=premnmx(XX);

YY=premnmx(YY);

YY

%創建網絡

net=newff(minmax(XX),[19,25,1],{'tansig','tansig','purelin'},'trainlm');

%下面使用遺傳算法對網絡進行優化

P=XX;

T=YY;

R=size(P,1);

S2=size(T,1);

S1=25;%隱含層節點數

S=R*S1+S1*S2+S1+S2;%遺傳算法編碼長度

aa=ones(S,1)*[-1,1];

popu=50;%種群規模

save data2 XX YY % 是將 xx,yy 二個變數的數值存入 data2 這個MAT-file，

initPpp=initializega(popu,aa,'gabpEval');%初始化種群

gen=100;%遺傳代數

%下面調用gaot工具箱，其中目標函數定義為gabpEval

[x,endPop,bPop,trace]=ga(aa,'gabpEval',[],initPpp,[1e-6 1 1],'maxGenTerm',gen,...

??'normGeomSelect',[0.09],['arithXover'],[2],'nonUnifMutation',[2 gen 3]);

%繪收斂曲線圖

figure(1)

plot(trace(:,1),1./trace(:,3),'r-');

hold on

plot(trace(:,1),1./trace(:,2),'b-');

xlabel('Generation');

ylabel('Sum-Squared Error');

figure(2)

plot(trace(:,1),trace(:,3),'r-');

hold on

plot(trace(:,1),trace(:,2),'b-');

xlabel('Generation');

ylabel('Fittness');

%下面將初步得到的權值矩陣賦給尚未開始訓練的BP網絡

[W1,B1,W2,B2,P,T,A1,A2,SE,val]=gadecod(x);

net.LW{2,1}=W1;

net.LW{3,2}=W2;

net.b{2,1}=B1;

net.b{3,1}=B2;

XX=P;

YY=T;

%設置訓練參數

net.trainParam.show=1;

net.trainParam.lr=1;

net.trainParam.epochs=50;

net.trainParam.goal=0.001;

%訓練網絡

net=train(net,XX,YY);


%程序二：適應值函數
function [sol, val] = gabpEval(sol,options)
% val - the fittness of this individual
% sol - the individual, returned to allow for Lamarckian evolution
% options - [current_generation]
load data2
nntwarn off
XX=premnmx(XX);
YY=premnmx(YY);
P=XX;
T=YY;
R=size(P,1);
S2=size(T,1);
S1=25;%隱含層節點數
S=R*S1+S1*S2+S1+S2;%遺傳算法編碼長度
for i=1:S,
? ?x(i)=sol(i);
end;
[W1, B1, W2, B2, P, T, A1, A2, SE, val]=gadecod(x);

%程序三：編解碼函數
function [W1, B1, W2, B2, P, T, A1, A2, SE, val]=gadecod(x)
load data2
nntwarn off
XX=premnmx(XX);
YY=premnmx(YY);
P=XX;
T=YY;
R=size(P,1);
S2=size(T,1);
S1=25;%隱含層節點數
S=R*S1+S1*S2+S1+S2;%遺傳算法編碼長度
% 前R*S1個編碼為W1
for i=1:S1,
? ? for k=1:R,
? ?? ?W1(i,k)=x(R*(i-1)+k);
? ? end
end
% 接著的S1*S2個編碼（即第R*S1個后的編碼）為W2
for i=1:S2,
? ?for k=1:S1,
? ?? ?W2(i,k)=x(S1*(i-1)+k+R*S1);
? ?end
end
% 接著的S1個編碼（即第R*S1+S1*S2個后的編碼）為B1
for i=1:S1,
? ?B1(i,1)=x((R*S1+S1*S2)+i);
end
% 接著的S2個編碼（即第R*S1+S1*S2+S1個后的編碼）為B2
for i=1:S2,
? ?B2(i,1)=x((R*S1+S1*S2+S1)+i);
end
% 計算S1與S2層的輸出
A1=tansig(W1*P,B1);
A2=purelin(W2*A1,B2);
% 計算誤差平方和
SE=sumsqr(T-A2);
val=1/SE; % 遺傳算法的適應值

針對第三問：該部分與問題二思路大致相同，但由于ADMET數據都是由0和1組成，因此若采用神經網絡模型時需要考慮隸屬度問題以及負值問題，同時根據建模常識盡量不要用同一種算法求解全部關鍵問題，但可以是對算法的升級，因此這里也可以采用深度學習相關的分類算法來做，主流的深度學習方法包括：AlexNet、VGG-Net、ResNet等；同時在這里推薦使用其它的分類算法，例如SVM支持向量機、隨機森林、adaboost、xgboost等。

對于第三問而言，主要是構建分類模型，這里有同學想的是構建五個不同的分類模型，比如綜合用神經網絡、SVM、、、等五個模型，其實完全沒必要，因為該分類問題目標類別只有兩種，但主流的分類方法都適合于多類別的劃分，因此在設計模型時（尤其是神經網絡模型，要注意隸屬度的問題），但在這里重點推薦邏輯回歸分類模型，目前Logistic回歸主要就是用在二分類領域。它是一個非線性的回歸模型， 其最大的好處恰恰是可以解決二元類問題，目前在很多行業，基本都是使用Logistic回歸來預判目標是好還是壞，因為它還彌補了其他黑盒模型(SVM、 神經網絡、隨機森林等)不具解釋性的缺點。

Logistic回歸分類代碼如下

分類案例：本文使用的例子是安德森鳶尾花卉數據集（iris），這個數據集中包含150個樣本，對應著它的150行，每一行包含這個樣本的4個特征（花萼長度、花萼寬度、花瓣長度、花瓣寬度）和樣本的類別標簽（0或1或2，它們分別代表不同的品種），所以數據集中存儲的是一個尺寸為150×5的二維矩陣（數據集中部分樣本的信息如下圖所示）。由于我們此次只研究二分類問題，所以我們只選用數據集的前100個樣本，即類別標簽為“0”或“1”的樣本。我們的任務是：建立一個邏輯回歸分類器，這個分類器可以通過樣本的4個特征值來預測樣本的標簽是“0”還是“1”。

?

iris數據集我已經上傳到百度云，有需要的小伙伴可以自行下載。

鏈接：https://pan.baidu.com/s/15p8-JOoui8taai0tSADNyA?
提取碼：k0cs?

以上數據集僅供參考，下面替換成自己的就行

% 讀取數據集，變量data存儲的是一個尺寸為150×5的矩陣

data = load('iris.data');

% 只取前100個樣本的信息，即前100行，變量useful_data存儲的是一個尺寸為100×5的矩陣

useful_data = data(1:100, :);

% 將特征與類別標簽分開存放

% 特征存放在變量X中，X存儲的是一個尺寸為100×4的矩陣

% 類別標簽存放在變量y中，y存儲的是一個尺寸為100×1的矩陣

X = useful_data(:, 1:4);

y = useful_data(:, 5);

% 變量m存儲的是變量X的行數，在這里為100

% 變量n存儲的是變量X的列數，在這里為4

[m,n] = size(X);

% 在變量X后加一列“1”，便于后面使用矩陣運算簡化步驟

X = [X ones(m, 1)];

% 初始化模型參數β=(ω,b)為0

beta = zeros(n+1, 1);

% 設置梯度下降迭代次數為1500次

iteration = 1500;

% 設置學習率為0.01

alpha = 0.01;

% 開始循環，用梯度下降更新參數

for iter = 1 : iteration

????z = X * beta;

????h = 1 ./ (1 + exp(-z));

????error = h - y;

????graident = X' * error;

????beta = beta - alpha / m * graident;

end

針對問題四：

問題4主要是探究哪些分子描述符對生物活性和ADMET效果更好，在這里我們可以按照兩步走的方法，首先是利用問題一選擇的20個變量和問題二構建的模型（在此以神經網絡模型為例），構建完神經網絡模型后，此時編寫相關的程序如下：

fun = @(x) -predict(net, x);

x0 = [0,0,0,0]; %你自己設置初值，（有可能需要轉置即x0 = [0,0,0,0]'）

[x, ymax] = fminsearch(fun, x0); x就是你要的位置

ymax = -ymax;%ymax就是極大值

此時就能找到能夠使生物活性達到最大值時上述分子描述符的取值，然后將上述初值帶入到問題三構建的分類模型中，看看ADMET分別取值是多少，若有三個性質較好則表示已經完成，若效果并不好，則返回第一步按照問題一中相關性高低原則逐漸剔除因子，指導既滿足分子活性最大又使得ADMET中至少三個性質較好。

?

總結

以上是生活随笔為你收集整理的2021年中国研究生数学建模竞赛D题参考思路的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。