2022年数模国赛C题（岭回归、区间预测、矩阵热力图、Fisher判别分类模型）——总结心得（附最后一次数模经历，Matlab\SPSS\Lingo的理解综合）_2022国赛c题第一问

一、国赛

第二次数模国赛
也就是昨天才结束的，选题很果断，选C，然后想半天，郭佬，给了个解释，很棒的解释，用SPSSPRO这个网站（大家可以去搜搜，很棒还有一个镝数图表，嗯数据给的少，放入Logistic逻辑回归的显著性都大于0.05，不构成影响，又代入什么逐步回归、岭回归模型，比较了一下，最终决定岭回归模型（矮子里的将军），最终每个化学元素都有对应的系数。预测方面，我们采用最值区间进行预测，在对于高钾和铅钡，二氧化硅、氧化钾，氧化钡都是未分化的区间分布较大，使用正态分布，选择最高的频数分布，作为预测区间

第二问
我们使用，Fisher判别分类模型，也是看的参考文献找到的模型（不得不说。我队长看文献找模型的能力真强），代入70%的训练集求出4个维度方程，

第三问
根据表三，求出质心坐标，只要落在质心附近就行，进行分类，说明范围即可（一定要查阅相关文献）

第四问
就是绘制散点图矩阵（超级棒的图），还有热力图，进行分析即可

想要总体思路图，可以私聊博主。

二、国赛代码展示

1.1 问题一

%% 读取
[num,txt,raw] = xlsread('附件.xlsx',1);
[num_g,txt_g,raw_g] = xlsread('高钾.xlsx');
[num_q,txt_q,raw_q] = xlsread('铅钡.xlsx');

%% 量化标识
for i = 2:59    
    for j = 2:5
        if cell2mat(txt(i,j)) == 'A'
           txt(i,j) = num2cell('1');
        elseif cell2mat(txt(i,j)) == 'B'
            txt(i,j) = num2cell('2');
        elseif cell2mat(txt(i,j)) == 'C'
            txt(i,j) = num2cell('3');
            
        elseif cell2mat(txt(i,3)) == '高钾'
            txt(i,j) = num2cell('1');
        elseif cell2mat(txt(i,3)) == '铅钡'
            txt(i,j) = num2cell('2');
        
        elseif isempty(cell2mat(txt(i,4)))
            txt(i,j) = num2cell('9');
        elseif cell2mat(txt(i,4)) == '黑'
            txt(i,j) = num2cell('1');
        elseif cell2mat(txt(i,4)) == '蓝绿'
            txt(i,j) = num2cell('2');
        elseif cell2mat(txt(i,4)) == '绿'
            txt(i,j) = num2cell('3');
        elseif cell2mat(txt(i,4)) == '浅蓝'
            txt(i,j) = num2cell('4');
        elseif cell2mat(txt(i,4)) == '浅绿'
            txt(i,j) = num2cell('5');
        elseif cell2mat(txt(i,4)) == '深蓝'
            txt(i,j) = num2cell('6');
        elseif cell2mat(txt(i,4)) == '深绿'
            txt(i,j) = num2cell('7');
        elseif cell2mat(txt(i,4)) == '紫'
            txt(i,j) = num2cell('8');
            
        elseif length(cell2mat(txt(i,5))) == 3 
            txt(i,j) = num2cell('1');
        elseif cell2mat(txt(i,5)) == '风化'
            txt(i,j) = num2cell('0');
        end
    end
end

%% 风化的数据
fenghua = [];
count = 0;
for i = 1:58
   if str2num(cell2mat(txt(i,5))) == 0
       count = count + 1;
       fenghua(count,1) = str2num(cell2mat(txt(i,2)));
       fenghua(count,2) = str2num(cell2mat(txt(i,3)));
       fenghua(count,3) = str2num(cell2mat(txt(i,4)));
   end
end

%% 风化画图
X = fenghua(:,3);  %颜色
Y = fenghua(:,2);  %类型
Z = fenghua(:,1);  %纹饰
% scatter3(X,Y,Z)
stem3(X,Y,Z)
set(gca,'xtick',0:1:9,'XTickLabel',{' ','黑','蓝绿','绿','浅蓝','浅绿','深蓝','深绿','紫','残缺'});
set(gca,'ytick',0:1:2,'YTickLabel',{' ','高钾','铅钡'});
set(gca,'ztick',0:1:3,'ZTickLabel',{' ','纹饰A','纹饰B','纹饰C'});

one_2_9 = 0;
one_2_4 = 0;
one_2_2 = 0;
one_2_1 = 0;
two_1_2 = 0;
three_2_2 = 0;
three_2_4 = 0;
three_2_5 = 0;
three_2_7 = 0;
three_2_8 = 0;
three_2_9 = 0;
for i = 1:33
    if fenghua(i,1) == 1
        if fenghua(i,3) == 9
            one_2_9 = one_2_9 + 1;
        elseif fenghua(i,3) == 4
            one_2_4 = one_2_4 + 1;
        elseif fenghua(i,3) == 2
            one_2_2 = one_2_2 + 1;
        elseif fenghua(i,3) == 1
            one_2_1 = one_2_1 + 1;
        end
    elseif fenghua(i,1) == 2
        two_1_2 = two_1_2 + 1;
    elseif fenghua(i,1) == 3
        if fenghua(i,3) == 2
            three_2_2 = three_2_2 + 1;
        elseif fenghua(i,3) == 4
            three_2_4 = three_2_4 + 1;
        elseif fenghua(i,3) == 5
            three_2_5 = three_2_5 + 1;
        elseif fenghua(i,3) == 7
            three_2_7 = three_2_7 + 1;
        elseif fenghua(i,3) == 8
            three_2_8 = three_2_8 + 1;
        elseif fenghua(i,3) == 9
            three_2_9 = three_2_9 + 1;
        end
    end
end
%% 风化的前提下，且玻璃为高钾
fenghua_gaojia = [];
count = 0;
for i = 1:33
   if fenghua(i,2) == 1
       count = count + 1;
       fenghua_gaojia(count,1) = fenghua(i,1);
       fenghua_gaojia(count,2) = fenghua(i,3);
   end
end

%% 风化的前提下，且玻璃为铅钡
fenghua_qianbei = [];
count = 0;
for i = 1:33
   if fenghua(i,2) == 2
       count = count + 1;
       fenghua_qianbei(count,1) = fenghua(i,1);
       fenghua_qianbei(count,2) = fenghua(i,3);
   end
end

%% 颜色与风化画图
yanse_fenghua = [2,0;9,6;0 1;12 8;1 2;0 2;4 3;2 2;4 0];
bar3(yanse_fenghua,1,'grouped')
set(gca,'ytick',1:1:9,'YTickLabel',{'黑','蓝绿','绿','浅蓝','浅绿','深蓝','深绿','紫','残缺'});
zlabel('计数','FontName','Song','FontSize', 14)
ylabel('颜色','FontName','Song','FontSize', 14)
legend('风化','无风化')

%% 纹饰与风化画图---叠
figure
wenshi_fenghua = [11 11;6 0;17 13];
bar3(wenshi_fenghua,0.4,'stacked')
set(gca,'ytick',1:1:3,'YTickLabel',{'纹饰A','纹饰B','纹饰C'});
zlabel('计数','FontName','Song','FontSize', 14)
ylabel('纹饰','FontName','Song','FontSize', 14)
legend('风化','无风化')

%% 类型与风化画图---叠
figure
leixing_fenghua = [6 12;28 12;];
bar3(leixing_fenghua,0.5,'stacked')
set(gca,'ytick',1:1:2,'YTickLabel',{'高钾','铅钡'});
zlabel('计数','FontName','Song','FontSize', 14)
ylabel('类型','FontName','Song','FontSize', 14)
legend('风化','无风化')

%% 纹饰与风化画图---并
figure
wenshi_fenghua = [11 11;6 0;17 13];
bar3(wenshi_fenghua,1,'grouped')
set(gca,'ytick',1:1:3,'YTickLabel',{'纹饰A','纹饰B','纹饰C'});
zlabel('计数','FontName','Song','FontSize', 14)
ylabel('纹饰','FontName','Song','FontSize', 14)
legend('风化','无风化')

%% 类型与风化画图---并
figure
leixing_fenghua = [6 12;28 12;];
bar3(leixing_fenghua,1,'grouped')
set(gca,'ytick',1:1:2,'YTickLabel',{'高钾','铅钡'});
zlabel('计数','FontName','Song','FontSize', 14)
ylabel('类型','FontName','Song','FontSize', 14)
legend('风化','无风化')

%% 高钾区间
g_no_winded = [];
g_yes_winded = [];
for i = 2:19
   if  cell2mat(raw_g(i,4)) == 0 % 未风化
       g_no_winded((length(g_no_winded) + 1)) = i;
   elseif  cell2mat(raw_g(i,4)) == 1 % 风化
       g_yes_winded((length(g_yes_winded) + 1)) = i;
   end
end

% 输入有无风化信息
g_no_winded_info = {};
g_yes_winded_info = {};

for i = 1:length(g_no_winded)
    g_no_winded_info(i,:) = raw_g(g_no_winded(i),:);
end

for i = 1:length(g_yes_winded)
    g_yes_winded_info(i,:) = raw_g(g_yes_winded(i),:);
end
% 最大最小
g_no_winded_max = [];
g_no_winded_min = [];

g_yes_winded_max = [];
g_yes_winded_min = [];
for j = 5:size(g_no_winded_info,2)
    g_no_winded_max(length(g_no_winded_max) + 1) = max(cell2mat(g_no_winded_info(:,j)));
    g_no_winded_min(length(g_no_winded_min) + 1) = min(cell2mat(g_no_winded_info(:,j)));
end

for j = 5:size(g_yes_winded_info,2)
    g_yes_winded_max(length(g_yes_winded_max) + 1) = max(cell2mat(g_yes_winded_info(:,j)));
    g_yes_winded_min(length(g_yes_winded_min) + 1) = min(cell2mat(g_yes_winded_info(:,j)));
end

% 获取最大减最小的范围
g_no_winded_range = g_no_winded_max - g_no_winded_min;
g_yes_winded_range = g_yes_winded_max - g_yes_winded_min;

g_range = g_no_winded_range - g_yes_winded_range

% 小画图展示 
X = 1:14;
subplot(2,1,1)
plot(X,g_no_winded_max);hold on
plot(X,g_no_winded_min);hold off

subplot(2,1,2)
plot(X,g_yes_winded_max);hold on
plot(X,g_yes_winded_min);hold off

%% 铅钡区间

q_no_winded = [];
q_yes_winded = [];
for i = 2:50
   if  cell2mat(raw_q(i,4)) == 0 % 未风化
       q_no_winded((length(q_no_winded) + 1)) = i;
   elseif  cell2mat(raw_q(i,4)) == 1 % 风化
       q_yes_winded((length(q_yes_winded) + 1)) = i;
   end
end

% 输入有无风化信息
q_no_winded_info = {};
q_yes_winded_info = {};

for i = 1:length(q_no_winded)
    q_no_winded_info(i,:) = raw_q(q_no_winded(i),:);
end

for i = 1:length(q_yes_winded)
    q_yes_winded_info(i,:) = raw_q(q_yes_winded(i),:);
end
% 最大最小
q_no_winded_max = [];
q_no_winded_min = [];

q_yes_winded_max = [];
q_yes_winded_min = [];
for j = 5:size(q_no_winded_info,2)
    q_no_winded_max(length(q_no_winded_max) + 1) = max(cell2mat(q_no_winded_info(:,j)));
    q_no_winded_min(length(q_no_winded_min) + 1) = min(cell2mat(q_no_winded_info(:,j)));
end

for j = 5:size(q_yes_winded_info,2)
    q_yes_winded_max(length(q_yes_winded_max) + 1) = max(cell2mat(q_yes_winded_info(:,j)));
    q_yes_winded_min(length(q_yes_winded_min) + 1) = min(cell2mat(q_yes_winded_info(:,j)));
end

% 获取最大减最小的范围
q_no_winded_range = q_no_winded_max - q_no_winded_min;
q_yes_winded_range = q_yes_winded_max - q_yes_winded_min;

q_range = q_no_winded_range - q_yes_winded_range

% 小画图展示 
X = 1:14;
subplot(4,1,1)
plot(X,g_no_winded_max);hold on
set(gca,'xtick',1:1:14,'XTickLabel',{'SiO2','Na2O','K2O','CaO','MgO','Al2O3','Fe2O3','CuO','PbO','BaO','P2O5','SrO','SnO2','SO2'});
plot(X,g_no_winded_min);hold off
set(gca,'xtick',1:1:14,'XTickLabel',{'SiO2','Na2O','K2O','CaO','MgO','Al2O3','Fe2O3','CuO','PbO','BaO','P2O5','SrO','SnO2','SO2'});

subplot(4,1,2)
plot(X,g_yes_winded_max);hold on
set(gca,'xtick',1:1:14,'XTickLabel',{'SiO2','Na2O','K2O','CaO','MgO','Al2O3','Fe2O3','CuO','PbO','BaO','P2O5','SrO','SnO2','SO2'});
plot(X,g_yes_winded_min);hold off
set(gca,'xtick',1:1:14,'XTickLabel',{'SiO2','Na2O','K2O','CaO','MgO','Al2O3','Fe2O3','CuO','PbO','BaO','P2O5','SrO','SnO2','SO2'});

subplot(4,1,3)
plot(X,q_no_winded_max);hold on
set(gca,'xtick',1:1:14,'XTickLabel',{'SiO2','Na2O','K2O','CaO','MgO','Al2O3','Fe2O3','CuO','PbO','BaO','P2O5','SrO','SnO2','SO2'});
plot(X,q_no_winded_min);hold off
set(gca,'xtick',1:1:14,'XTickLabel',{'SiO2','Na2O','K2O','CaO','MgO','Al2O3','Fe2O3','CuO','PbO','BaO','P2O5','SrO','SnO2','SO2'});

subplot(4,1,4)
plot(X,q_yes_winded_max);hold on
set(gca,'xtick',1:1:14,'XTickLabel',{'SiO2','Na2O','K2O','CaO','MgO','Al2O3','Fe2O3','CuO','PbO','BaO','P2O5','SrO','SnO2','SO2'});
plot(X,q_yes_winded_min);hold off
set(gca,'xtick',1:1:14,'XTickLabel',{'SiO2','Na2O','K2O','CaO','MgO','Al2O3','Fe2O3','CuO','PbO','BaO','P2O5','SrO','SnO2','SO2'});
%% 铅钡_风化/未风化__二氧化硅和氧化铅
clc

q_no_winded_SiO2 = cell2mat(q_no_winded_info(:,5));
q_yes_winded_SiO2 = cell2mat(q_yes_winded_info(:,5));
q_no_winded_PbO = cell2mat(q_no_winded_info(:,13));
q_yes_winded_PbO = cell2mat(q_yes_winded_info(:,13));

g_no_winded_SiO2 = cell2mat(g_no_winded_info(:,5));
g_yes_winded_SiO2 = cell2mat(g_yes_winded_info(:,5));

% 画频率分布直方图
[counts,centers] = hist(g_yes_winded_SiO2, 5);
figure
bar(centers, counts / sum(counts))

% 分布参数拟合
[mu,sigma]=normfit(g_yes_winded_SiO2);

% 画已知分布的概率密度曲线
x1 = 88:0.1:100;
y1 = pdf('Normal', x1, mu,sigma);
hold on
plot(x1, y1*2.5)
title('g-no-winded-SiO_2')



%% 预测与真实
X = 1:18;
Y1 = [0 0 0 0 0 0 0 1 1 1 1 0 0 0 0 0 1 1];
Y2 = [0.074 0.494 -0.097 -0.029 -0.011 0.023 -0.069 0.881 0.914 0.965 0.904 0.028 -0.116 0.081 0.055 0.239 0.788 0.875];
plot(X,Y1); hold on;
plot(X,Y2); hold off;

%% 岭回归
x1 = cell2mat(raw_q(2:end,5));x2 = cell2mat(raw_q(2:end,6));
x3 = cell2mat(raw_q(2:end,7));x4 = cell2mat(raw_q(2:end,8));
x5 = cell2mat(raw_q(2:end,9));x6 = cell2mat(raw_q(2:end,10));
x7 = cell2mat(raw_q(2:end,11));x8 = cell2mat(raw_q(2:end,12));
x9 = cell2mat(raw_q(2:end,13));x10 = cell2mat(raw_q(2:end,14));
x11 = cell2mat(raw_q(2:end,15));x12 = cell2mat(raw_q(2:end,16));
x13 = cell2mat(raw_q(2:end,17));x14 = cell2mat(raw_q(2:end,18));

y = cell2mat(raw_q(2:end,4));


x = [x1,x2,x3,x4,x5,x6,x6,x7,x8,x9,x10,x11,x12,x13,x14];
% 标准化
x = (x - min(x))./(max(x) - min(x));
%先画出岭迹图，以便选取合适的岭参数
k = 0:1e-3:10; %岭参数
b = ridge(y,x,k); %回归系数
figure
plot(k,b,'LineWidth',1.5)
set(gca,'xtick',10:10:100,'XTickLabel',{'0.1','0.2','0.3','0.4','0.5','0.6','0.7','0.8','0.9','1'});
xlabel('k')
ylabel('β')

%每个k对应的残差平方和的图，要选取会使残差平方和小的k值
knum=size(b,2);
sse=zeros(knum);
y_gu=zeros(49,knum);
for j=1:knum
    t=0;
for i=1:49
    y_gu(i,j)=b(1,j)*x1(i)+b(2,j)*x2(i)+b(3,j)*x3(i)+b(4,j)*x4(i)+b(5,j)*x5(i)+b(6,j)*x6(i)+b(7,j)*x7(i);
    t=(y_gu(i,j)-y(i))^2+t;
end
    sse(j)=t;
end

figure
plot(k,sse(:,1),'LineWidth',2)
xlabel('k')
ylabel('SSE')
title('残差平方和与k的关系图')
% exportgraphics(gcf,['C:\Users\DELL\Desktop\图\figure',num2str(i),'.jpg'])
warning off


%% 高钾的平均值
g_fenghua = [];
g_weifenghua = [];

Y_fenghua = mean(g_fenghua);
Y_weifenghua = mean(g_weifenghua);

X = 1:7;
plot(X,Y_fenghua);hold on
plot(X,Y_weifenghua);hold off

legend('风化','未分化')




1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396

1.2 问题二

%% 读取数据 '高钾_铅钡、风化_未风化.xlsx'
[num_g_yes_winded,txt_g_yes_winded,raw_g_yes_winded] = xlsread('高钾_铅钡、风化_未风化.xlsx',1);
[num_g_no_winded,txt_g_no_winded,raw_g_no_winded] = xlsread('高钾_铅钡、风化_未风化.xlsx',2);
[num_q_yes_winded,txt_q_yes_winded,raw_q_yes_winded] = xlsread('高钾_铅钡、风化_未风化.xlsx',3);
[num_q_no_winded,txt_q_no_winded,raw_q_no_winded] = xlsread('高钾_铅钡、风化_未风化.xlsx',4);

[num_2,txt_2,raw_2] = xlsread('附件.xlsx',2);
%% 求均值
g_yes_winded_mean = mean(cell2mat(raw_g_yes_winded(2:end,5:end)));
g_no_winded_mean = mean(cell2mat(raw_g_no_winded(2:end,5:end)));
q_yes_winded_mean = mean(cell2mat(raw_q_yes_winded(2:end,5:end)));
q_no_winded_mean = mean(cell2mat(raw_q_no_winded(2:end,5:end)));

all_mean = [g_yes_winded_mean;g_no_winded_mean;q_yes_winded_mean;q_no_winded_mean]';
all_mean = [all_mean(1:end-2,:)] % 删除低于20%的
all_mean = [all_mean(1,:);all_mean(3:end,:)]; % 删除25%的


all_mean = all_mean'
%% Fisher判别法
all_mean = all_mean'
x1=all_mean(1:5,:);
x2=all_mean(6:9,:);
sample=all_mean(10:11,:);
y=fisher(x1,x2,sample)

%% 求坐标 求距离 求准确度
% 获取需要数据并删除15和17行和氧化钠和最后两列
X = [raw_2(2:18,2),raw_2(2:18,4:end-3);
     raw_2(20,2),raw_2(20,4:end-3);
     raw_2(22:end,2),raw_2(22:end,4:end-3)];
X = cell2mat(X);
% 清除NaN
for i = 1:size(X,1)
   for j = 1:size(X,2)
      if isnan(X(i,j))
          X(i,j) = 0;
      end
   end
end

% 选取50组往后的数据
X = X(51:end,:)

for i = 1:size(X,1)
   F(i,1) = 0.065*X(i,1)+0.253*X(i,2)-0.211*X(i,3)-0.043*X(i,4)-17.985*X(i,5)+0.031*X(i,6)-1.518*X(i,7)-0.031*X(i,8)+0.141*X(i,9)+0.002*X(i,10)+4.491*X(i,11)-3.62;
   F(i,2) = -0.045*X(i,1)+0.541*X(i,2)-0.324*X(i,3)-0.043*X(i,4)-6.559*X(i,5)-0.102*X(i,6)+1.511*X(i,7)+0.013*X(i,8)+0.017*X(i,9)+0.459*X(i,10)+2.195*X(i,11)-4.14;
   F(i,3) = 0.102*X(i,1)+0.206*X(i,2)-0.505*X(i,3)-0.271*X(i,4)+13.49*X(i,5)-0.061*X(i,6)+3.054*X(i,7)+0.012*X(i,8)+0.115*X(i,9)+0.367*X(i,10)-0.511*X(i,11)-4.343;
   F(i,4) = 0.053*X(i,1)-0.015*X(i,2)+0.115*X(i,3)+0.18*X(i,4)-1.645*X(i,5)-0.091*X(i,6)-1.403*X(i,7)-0.166*X(i,8)-0.056*X(i,9)+0.399*X(i,10)+4.631*X(i,11)-0.663;
end





1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55

1.3 问题三

%% 读取数据
[num,txt,raw] = xlsread('附件.xlsx',3);

%% 计算坐标
% 数据清洗
X = [cell2mat(raw(2:end,3)),cell2mat(raw(2:end,5:end-2))];
% NaN变0
for i = 1:size(X,1)
   for j = 1:size(X,2)
      if isnan(X(i,j))
          X(i,j) = 0;
      end
   end
end

for i = 1:size(X,1)
   F(i,1) = 0.065*X(i,1)+0.253*X(i,2)-0.211*X(i,3)-0.043*X(i,4)-17.985*X(i,5)+0.031*X(i,6)-1.518*X(i,7)-0.031*X(i,8)+0.141*X(i,9)+0.002*X(i,10)+4.491*X(i,11)-3.62;
   F(i,2) = -0.045*X(i,1)+0.541*X(i,2)-0.324*X(i,3)-0.043*X(i,4)-6.559*X(i,5)-0.102*X(i,6)+1.511*X(i,7)+0.013*X(i,8)+0.017*X(i,9)+0.459*X(i,10)+2.195*X(i,11)-4.14;
   F(i,3) = 0.102*X(i,1)+0.206*X(i,2)-0.505*X(i,3)-0.271*X(i,4)+13.49*X(i,5)-0.061*X(i,6)+3.054*X(i,7)+0.012*X(i,8)+0.115*X(i,9)+0.367*X(i,10)-0.511*X(i,11)-4.343;
   F(i,4) = 0.053*X(i,1)-0.015*X(i,2)+0.115*X(i,3)+0.18*X(i,4)-1.645*X(i,5)-0.091*X(i,6)-1.403*X(i,7)-0.166*X(i,8)-0.056*X(i,9)+0.399*X(i,10)+4.631*X(i,11)-0.663;
end

%% 绘制距离图
scatter(-1.4,-1.7,150,'filled','d');hold on %高未
scatter(-0.8,1.2,150,'filled','h');hold on %高风
scatter(2.6,0.3,150,'filled','s');hold on  %铅未
scatter(3.5,3.7,150,'filled','p');hold on  %铅风

scatter(-1.2,-1.5,100,'d');hold on
X = [-0.4,-0.6];
Y = [1.4,2.0];
scatter(X,Y,100,'h','r');hold on
% scatter(-0.4,1.4,40,'h','r');hold on
% scatter(-0.6,2.0,40,'h','r');hold on

X = [2.3 2.6 3.0];
Y = [0.6 1.1 -0.1];
scatter(X,Y,100,'s','k');hold on
% scatter(2.3,0.6,40,'filled','s','g');hold on
% scatter(2.6,1.1,40,'filled','s','g');hold on
% scatter(3.0,-0.1,40,'filled','s','g');hold on


X = [3.3 2.5];
Y = [2.7 3.7];
scatter(X,Y,100,'p');hold on
xlim([-2,4]);
ylim([-2,4]);
lgd = legend('高钾未风化质心','高钾风化质心','铅钡未风化质心','铅钡风化质心')

%% 45样本

scatter(-1.6,-1.2,150,'filled','d');hold on %高未
scatter(-0.82,1.4,150,'filled','h');hold on %高风
scatter(2.7,0.5,150,'filled','s');hold on  %铅未
scatter(3.2,3.3,150,'filled','p');hold on  %铅风

scatter(-1.2,-1.5,100,'d');hold on
X = [-0.4,-0.6];
Y = [1.4,2.0];
scatter(X,Y,100,'h','r');hold on
% scatter(-0.4,1.4,40,'h','r');hold on
% scatter(-0.6,2.0,40,'h','r');hold on

X = [2.3 2.6 3.0];
Y = [0.6 1.1 -0.1];
scatter(X,Y,100,'s','k');hold on
% scatter(2.3,0.6,40,'filled','s','g');hold on
% scatter(2.6,1.1,40,'filled','s','g');hold on
% scatter(3.0,-0.1,40,'filled','s','g');hold on


X = [3.3 2.5];
Y = [2.7 3.7];
scatter(X,Y,100,'p');hold on
xlim([-2,4]);
ylim([-2,4]);
lgd = legend('高钾未风化质心','高钾风化质心','铅钡未风化质心','铅钡风化质心')
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78

1.4 问题四

A = []
A = [A(:,1),A(:,4),A(:,6:end-3)];
[S,AX,BigAx,H,HAx] = plotmatrix(A)


%%
[num,txt,raw] = xlsread('gaojia.xlsx');

for i = 2:13
   for j = 2:13 
       tmp = strsplit(cell2mat(raw(i,j)),'(');
       tmp = tmp(1);
       ch = cell2mat(tmp);
       doub = str2num(ch);
       A(i-1,j-1) = doub;
   end
end

xvalues = {'SiO_2', 'Na_2O', 'CuO', 'P_2O_5', 'K_2O','CaO','MgO','Fe_2O_3','PbO','Al_2O_3','SrO','BaO'};
yvalues = fliplr(xvalues)
h = heatmap(xvalues,xvalues,A);



1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25

三、对于软件的理解

3.1 Matlab

一定要多去看看函数，为什么野路子和科班出身的差距在哪？野路子是有个函数你甚至都不知道有这个函数可以用，科班出身就是函数都差不多知道了，要用的时候去查查就知道了，这一点真的很恐怖

3.1.1 表格的读取

[num,txt,raw] = xlsread(‘附件.xlsx’,1);
这个函数超级好用，读表格就选它，raw是最好的，毕竟他把所有的数据形式，都变成元组形式，其实一开始，我也就得不咋地好用，不过呢，用多了，就知道它的香了给出相关函数

3.1.2 元胞数组的相关函数的转换

1、cell2mat：将cell转换为mat的char型
2、str2num：将mat从char转换为double型
3、cellstr：将char转cell
4、num2str：将double转char
5、num2cell：将double直接转cell

对于第一个cell2mat嗯，不一定是cell转char，如果里面是数组，那直接就是数组了，如果元胞里面是数字的话，你想调出这个数字，步骤如下，先cell2mat变成char型，再str2num转成double形式就行，恩方便

3.1.3 图像方面

figure建立图窗

figure(1)
x = 1:10:100;
y = [20 30 45 40 60 65 80 75 95 90];
yneg = [1 3 5 3 0 0 0 0 3 3];
ypos = [2 5 3 5 2 5 2 2 5 5];
xneg = [1 3 5 3 5 3 6 4 3 3];
xpos = [2 5 3 5 2 5 2 2 5 5];

h = errorbar(x,y,yneg,ypos,xneg,xpos,'LineStyle','--','Marker','s','MarkerSize',5,'MarkerEdgeColor','red','MarkerFaceColor','red');
set(gca,'xtick',0:25:100,'XTickLabel',{'00:00','06:00','12:00','18:00','24:00'});
set(gca,'FontName','times new roman','FontSize', 14);
xlabel('你好啊','FontName','Song','FontSize', 14)
ylabel('你好啊','FontName','Song','FontSize', 14)
box off
title('\it\fontname{Times New Roman}x\rm\fontname{Times New Roman}\zeta')
xlabel('\it\fontname{Times New Roman}x','FontSize',20)
ylabel('\it\fontname{Times New Roman}figure\sl(\it\fontname{Times New Roman}x\rm)')
set(h,'LineWidth',2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18

或者在查看菜单栏下，有个属性编辑器，你自己点点就知道干啥了


复制图窗可以最大程度上保留图像像素

介绍一个基础画图比较全面的系列文章MATLAB科学绘图

3.1.4 判空

isnan和isempty，嗯，用过的都说好用，你们可以查查怎么用，做数据缺失时，贼好用

3.1.5 对于有字符串影响的处理

比方说我队长让我出个热力图，数据都出来了，

我一看傻眼了，后面圆括号里是啥啊，嗯，需要去了，总不可能一个个手输入，

[num,txt,raw] = xlsread('gaojia.xlsx');

for i = 2:13
   for j = 2:13 
       tmp = strsplit(cell2mat(raw(i,j)),'(');
       tmp = tmp(1);
       ch = cell2mat(tmp);
       doub = str2num(ch);
       A(i-1,j-1) = doub;
   end
end
1
2
3
4
5
6
7
8
9
10
11

strsplit这个小函数还是特别好用的

3.1.6 矩阵翻转和旋转

Matlab对矩阵进行旋转、左右上下翻转、重新排列

3.2 Lingo

就是看哔站啦，看看视频，去学习，不过我记得几个注意事项吧，Lingo里面的灵敏度检验不要用，因为完全不知道里面内部是怎么实现的，第二，就是<=和>=好像在Lingo里面全部都是单纯的大于和小于

3.3 SPSS

SPSS的话，自己多去点点，或者像之前介绍的SPSSPRO一样也可以,点点

3.4 Eviews

时间序列的神器，嗯，不是很会用，B站上都有

3.5 亿图

论文的yyds，反正不是很会用，哈哈哈哈

3.6 Latex

emmm，论文排版神器，结果还是没去学

四、给初学者一些建议（针对软件方向）

前期你可以用SPSS点点，Lingo写写优化，Eviews做时间序列，但是你必须要会用Matlab，多多去尝试Matlab，当你会用Matlab的时候，其实其他的软件就不用再去深究了。Matlab的话，可以多去看看CSDN上的，很多轮子都造好了，比如说这次Fisher判别式，CSDN的代码都是一大堆的。emmm，遇到不会的，一定先自己在网上找办法，不会了，再去问学长学姐，或者同级之间的那些大牛，去问他们，他们的见识和理解一定可以帮到你的，还有一点就是。如果遇到自己不会的函数，去看官网，官网一定是最好最权威的解释，一定更全面。对于软件，模型的知识一定要懂，模型检验的话，灵敏度检验yyds，就是外部检验对整个模型的检验。
出的图，一定要去打磨！打磨！！一张图是否好看，就是看它信息量大不大，信息量越大，图就越好看，没人会觉得一个简单的折线图有多喊看的（这个点，我也是最近才从一个模型大佬那知道的，惭愧了，不过真的很有效）。就这些吧，嗯，一定要配合好模型队长，记住，你可以不会，但你一定要去学，9点钟回宿舍了，去学怎么写，我记得有一次模拟赛，我Lingo不会写，嗯，选题的当晚，回到宿舍看是逛哔站，一直学，学到2点，第二天就很有信心站在队长面前说，我可以解出来。
软件最厉害的地方在于，模型你尽管建，软件我都能求出来，哈哈哈哈，加油了，各位小软件~~~

附、一些话吧

作为一个软件，和模型队长也比赛很长时间了呀，从大一一直到大三，如果算上这次数模国赛的话，总共有6次比赛了，加上模拟赛就是14次数模论文写作，通了7次宵，两次五一赛，两次数模国赛，一次亚太赛，一次美赛，还有一次数模的国赛模拟总共7次通宵。昨天比完国赛，队长可能不再参赛了。很难受吧，毕竟快两年的战友了。我也一直很幸运遇到我的队长，每次都能悬悬的完赛，也不容易，哈哈哈哈。论文换了一次。也听说是队长的舍友选队友的时候，是跟我一个家乡的，就选了我。而我当时选队长时，嗯，4–6个模型找我，我就在想，我一定要找个城南校区（我是中心校区的），因为感觉中心模型太菜了，所以想找个城南的，之所以我选软件，原因有二，其一，是因为我这个专业需要学Matlab（最主要），其二，我是班长，不想再负责什么了（当时的我真的特别累吧），因为我们是选了模型就相当于当队长。
下面，我尽量就按照时间线回忆吧

第一次小模拟赛（好像是校赛）
算是协会对新队伍，第一个磨合测试题，嗯，我记得当时我用Lingo解出一个规划方程，好像是最短路径问题，我记得我当时特别激动反正，就是害怕自己不够用，哈哈哈哈，其余的，我记不太清了

第一次，五一赛
当时，我第一次见我队长，挺漂亮？还行吧，瘦瘦的，高高的，第一次，emmm，哈哈哈哈，我记得我当时贼怂，因为软件嘛，不咋地，还是需要发育的，怂怂的，反正我特怂（我感觉我队长肯定对我不是这么想的），选题啥A题，疫苗那题，优化，嗯，差点搞死我（哭了），我真的算不出来那个啥啥的时间序列。有个直系学长，在我后面，已经将全部代码都写出来了，正在跑代码。好吧，我当时觉得自己菜的要死，回去都不好意思见队长，哈哈哈哈，比赛的每个晚上，我应该都没在休息，都是想尽自己的力，把自己的思路啊，代码啊，写下来，一直都是这样，因为不好意思偷懒，哈哈哈哈。人生第一次通宵，也就是在那时，出现了生理反应，想干呕，不适等，最终也完赛了

暑假留校的八次模拟赛
暑假留校了啊，挺开心，遇到了很好的舍友，一个大二，两个大三，三个大一，大一的就是我，还有两个是城南的，一个’大佐’，'一个茂盛哥，大三是来考研的，其中一个中途走了，家里有事，另一个是应物学长，大二的那个学长，机电的，人也贼好，反正舍友都贼好。嗯，学长都特别好，我记得刚开始吧，我们去二期机房比赛的，哈哈哈哈，分配给我们队一个大二队伍来指导我们，一个软件学姐来指导我（嗯，跟那个学姐在暑假快结束了，白了一下，谈了一年多，分了），认识了很多人，董佬，徐佬，吕佬（有两个一个中心，一个城南），薛佬。女生的话，emmm认识，但都没有过多的交集，也就是第一次模拟赛吧，论文退了，身体不舒服了，换了一个，并且也一直跟我们比赛到了现在。早起跑步，打篮球，暴风雨来的时候，扫荡地超，1号楼宿舍旁的树倒了那晚。写了8次模拟赛，很累，最终我们组差不多就是物理题直接不看，优化题可以考虑，评价类首选。也想起来，徐州疫情严重，我们都在第一食堂三楼，聚在一起写题，嗯，早上六点的太阳，我那一个月真的是天天见啊，占个座，就去晨跑了，队长嘛，人很好，很负责，很幸运。在这里，我真的很感谢一个人，董佬，他总是能认真的回答我的每一个问题，真的，他人很好，很赞。哈哈哈哈，数模的每一个人都很棒，顶呱呱的。真的是很棒的一次经历，数模协会都是很优秀的存在。对了，忘说了一点，我是个急性子，我队长是慢性子，不过还好，最后总是能突破自己，写完题目，贼强。我也每次都挖苦队长，针对我，哈哈哈哈，我也每次只是开个玩笑，希望队长不要放心上，也经常和队长吵架，题目的理解问题情况，还有和队长的磨合出现的摩擦吧

亚太赛和美赛
都是线上的比的，记不太清，反正美赛我记得是和祖老一起在休息室通宵的

第二次五一赛
emmm，我也记不太清了，继续线上通宵,哭，这次和祖老与吕佬一起过的

总的来说，挺不错，很棒的数模协会，很棒的数模人，很棒的自己，很棒的经历