因子分析

# Conflicts:
#	盘荣博/数据可视化/week1.py

.idea/.gitignore

@@ -0,0 +1,3 @@
+# 默认忽略的文件
+/shelf/
+/workspace.xml

.idea/inspectionProfiles/Project_Default.xml

@@ -0,0 +1,12 @@
+<component name="InspectionProjectProfileManager">
+  <profile version="1.0">
+    <option name="myName" value="Project Default" />
+    <inspection_tool class="PyPackageRequirementsInspection" enabled="true" level="INFORMATION" enabled_by_default="true">
+      <option name="ignoredPackages">
+        <value>
+          <list size="0" />
+        </value>
+      </option>
+    </inspection_tool>
+  </profile>
+</component>

.idea/inspectionProfiles/profiles_settings.xml

@@ -0,0 +1,6 @@
+<component name="InspectionProjectProfileManager">
+  <settings>
+    <option name="USE_PROJECT_PROFILE" value="false" />
+    <version value="1.0" />
+  </settings>
+</component>

.idea/misc.xml

@@ -0,0 +1,4 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<project version="4">
+  <component name="ProjectRootManager" version="2" project-jdk-name="Python 3.11 (mycode)" project-jdk-type="Python SDK" />
+</project>

.idea/modules.xml

@@ -0,0 +1,8 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<project version="4">
+  <component name="ProjectModuleManager">
+    <modules>
+      <module fileurl="file://$PROJECT_DIR$/.idea/mycode.iml" filepath="$PROJECT_DIR$/.idea/mycode.iml" />
+    </modules>
+  </component>
+</project>

.idea/mycode.iml

@@ -0,0 +1,8 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<module type="PYTHON_MODULE" version="4">
+  <component name="NewModuleRootManager">
+    <content url="file://$MODULE_DIR$" />
+    <orderEntry type="inheritedJdk" />
+    <orderEntry type="sourceFolder" forTests="false" />
+  </component>
+</module>

.idea/vcs.xml

@@ -0,0 +1,6 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<project version="4">
+  <component name="VcsDirectoryMappings">
+    <mapping directory="" vcs="Git" />
+  </component>
+</project>

ACM/DongDong认亲戚.py

@@ -0,0 +1,2 @@
+a=[0]*int(2e4)
+print(len(a))

ACM/World Finals.py

@@ -0,0 +1,54 @@
+from functools import cmp_to_key
+c=[[],[]]
+string=[[],[]]
+for i in range(2):
+    a = eval(input())
+    while a>0:
+        s=input().split()
+        s[1],s[2]=eval(s[1]),eval(s[2])
+        c[i].append(s)
+        string[i].append(s[0])
+        a-=1
+def cmp(a,b):
+    if a[1]>b[1]:
+        return -1
+    elif a[1]==b[1]:
+        if a[2]<b[2] :
+            return -1
+        else:
+            return 1
+    else :
+        return 1
+c[0]=sorted(c[0],key=cmp_to_key(cmp))
+# print(c[0])
+c[1]=sorted(c[1],key=cmp_to_key(cmp))
+# print(c[1])
+# m='123'
+# n='123'
+# print(m==n)
+st="lzr010506"
+s=set(string[0])&set(string[1])
+# for i in c[0]:
+#     for j in c[1]:
+#         if i[0] == j[0] :
+#             s.append(j[0])
+# print(s)
+ans1=0
+ans2=0
+i=0
+while c[0][i][0] != st:
+    # print(c[0][i][0])
+    if c[0][i][0] in s:
+        ans1+=1
+    i+=1
+ans1= i+1-ans1
+i=0
+while c[1][i][0]!=st:
+    # print(c[1][i][0])
+    if c[1][i][0] in s:
+        ans2+=1
+    i+=1
+ans2=i+1-ans2
+# print(ans1,ans2)
+print(min(ans1,ans2))
+# print(c[0],'\n',c[1])

ACM/[NOIP2004]合并果子.py

@@ -0,0 +1,18 @@
+import heapq
+import ctypes
+a=int(input())
+b=input().split()
+for i in range(len(b)):
+    b[i]=int(b[i])
+heapq.heapify(b)
+# print(b)
+sum=0
+while len(b) != 1:
+    c,d=heapq.nsmallest(2,b)
+    # heapq.heappop(b)
+    # heapq.heappop(b)
+    b.remove(c)
+    b.remove(d)
+    sum+=c+d
+    heapq.heappush(b,c+d)
+print(sum)

ACM/值周.py

@@ -0,0 +1,28 @@
+import functools
+a,b =map(int, input().split(' '))
+nums=[]
+for i in range(b):
+    m,n=map(int,input().split(' '))
+    if m>n:continue;
+    nums.append([m,n])
+def cmp(a,b):#升序
+    if(a[0]>b[0]):
+        return 1
+    else :
+        return -1#不变
+nums=sorted(nums,key=functools.cmp_to_key(cmp))
+# for i in nums:
+#     print(i)
+sum = 0
+start,end=nums[0][0],nums[0][1]
+count=1
+while count<len(nums):
+    if nums[count][0]>=start and nums[count][0]<=end:
+        if end<=nums[count][1]:
+            end=nums[count][1]
+    else :
+        sum+=end-start+1
+        start=nums[count][0];end=nums[count][1]
+    count+=1
+sum+=end-start+1
+print(a-sum+1)

ACM/校门外的树.py

@@ -0,0 +1,7 @@
+a,b =map(int, input().split(' '))
+nums=[]
+for i in range(b):
+    m,n=map(int,input().split(' '))
+    nums.extend(range(m,n+1))
+nums=set(nums)
+print(a+1-len(nums))

姚安欣/AHP.m

@@ -0,0 +1,64 @@
+%AHP步骤
+
+clc,clear,close all;
+A=[1,2,3,5
+   1/2,1,1/2,2
+   1/3,2,1,2
+   1/5,1/2,1/2,1];
+[row,col]=size(A);
+
+%判断矩阵一致性检验
+n=col;
+maxlam=max(eig(A));
+RI=[0 0 0.58 0.90 1.12 1.24 1.32 1.41 1.45];
+CI=(maxlam-n)/(n-1);
+CR=CI/RI(n);
+
+%判断矩阵确定权重
+for i=1:col
+    sumcol=sum(A(:,i));
+    for j=1:row
+        A(j,i)=A(j,i)/sumcol;
+    end
+end
+
+weig=zeros(row,1);
+for i=1:row
+    sumrow=sum(A(i,:));
+    weig(i)=sumrow/n;
+end
+
+%各个指标归一化   按列单位化
+data= [1686.4 3183 12000 397
+       903.6 1916.4 3439.6 43
+       837.6 817.6 4748 1159
+       824.9 1296.4 12000 442
+       2110.2 1465.7 6199.5 228]; 
+[rowd,cold]=size(data);
+for i=1:cold
+    sumcold=sum(data(:,i));
+    for j=1:rowd
+        data(j,i)=data(j,i)/sumcold;
+    end
+end
+
+
+%按权重计算分数
+
+score=data*weig;
+
+
+projectNames={'老番茄','何同学','木鱼水心','凉风','罗翔'};
+figure; 
+bar(score);%条形图
+
+set(gca, 'XTickLabel', projectNames);  %每个条形图标签
+xlabel('博主');  
+ylabel('加权总分'); 
+title('得分');  
+  
+grid on; % 网格线
+
+
+
+

姚安欣/LPej1_3.m

@@ -1,40 +1,40 @@
-clc; clear; close all; format long g;
+clc; clear; close all; format long g;
-prob=optimproblem;
+prob=optimproblem;
-x=optimvar('x',9,1,'LowerBound',0,'Type','integer');
+x=optimvar('x',9,1,'LowerBound',0,'Type','integer');
-p=[0.25,0.35,0.50];
+p=[0.25,0.35,0.50];
-r=[1.25,2.00,2.80];
+r=[1.25,2.00,2.80];
-ot=[5,10,0;7,9,12;6,8,0;4,0,11;7,0,0];
+ot=[5,10,0;7,9,12;6,8,0;4,0,11;7,0,0];
-total=1000:200:3000;
+total=1000:200:3000;
-odp=[300/6000,321/10000,250/4000,783/7000,200/4000];
+odp=[300/6000,321/10000,250/4000,783/7000,200/4000];
-
+
-X=[];
+X=[];
-Q=[];
+Q=[];
-
+
-for i=1:length(total)
+for i=1:length(total)
-   % prob.Objective=(ot(1,1)*ones(1,3)*x(1:3)+ot(1,2)*x(7))*odp(1)+(ot(2,1)*ones(1,3)*x(4:6)+ot(2,2)*x(8)+ot(2,3)*x(9))*odp(2)+(ot(3,1)*(x(1)+x(4))+ot(3,2)*(x(7)+x(8)))*odp(3)+(ot(4,1)*(x(2)+x(5))+ot(4,3)*x(9))*odp(4)+(ot(5,1)*(x(3)+x(6)))*odp(5)-(r(1)-p(1))*ones(1,6)*x(1:6)-(r(2)-p(2))*ones(1,2)*x(7:8)-(r(3)-p(3))*x(9);
+   % prob.Objective=(ot(1,1)*ones(1,3)*x(1:3)+ot(1,2)*x(7))*odp(1)+(ot(2,1)*ones(1,3)*x(4:6)+ot(2,2)*x(8)+ot(2,3)*x(9))*odp(2)+(ot(3,1)*(x(1)+x(4))+ot(3,2)*(x(7)+x(8)))*odp(3)+(ot(4,1)*(x(2)+x(5))+ot(4,3)*x(9))*odp(4)+(ot(5,1)*(x(3)+x(6)))*odp(5)-(r(1)-p(1))*ones(1,6)*x(1:6)-(r(2)-p(2))*ones(1,2)*x(7:8)-(r(3)-p(3))*x(9);
-    dp1=(ot(1,1)*ones(1,3)*x(1:3)+ot(1,2)*x(7))*odp(1);
+    dp1=(ot(1,1)*ones(1,3)*x(1:3)+ot(1,2)*x(7))*odp(1);
-    dp2=(ot(2,1)*ones(1,3)*x(4:6)+ot(2,2)*x(8)+ot(2,3)*x(9))*odp(2);
+    dp2=(ot(2,1)*ones(1,3)*x(4:6)+ot(2,2)*x(8)+ot(2,3)*x(9))*odp(2);
-    dp3=(ot(3,1)*(x(1)+x(4))+ot(3,2)*(x(7)+x(8)))*odp(3);
+    dp3=(ot(3,1)*(x(1)+x(4))+ot(3,2)*(x(7)+x(8)))*odp(3);
-    dp4=(ot(4,1)*(x(2)+x(5))+ot(4,3)*x(9))*odp(4);
+    dp4=(ot(4,1)*(x(2)+x(5))+ot(4,3)*x(9))*odp(4);
-    dp5=(ot(5,1)*(x(3)+x(6)))*odp(5);
+    dp5=(ot(5,1)*(x(3)+x(6)))*odp(5);
-    dp=dp1+dp2+dp3+dp4+dp5;
+    dp=dp1+dp2+dp3+dp4+dp5;
-    prob.Objective=dp-(r(1)-p(1))*ones(1,6)*x(1:6)-(r(2)-p(2))*ones(1,2)*x(7:8)-(r(3)-p(3))*x(9);
+    prob.Objective=dp-(r(1)-p(1))*ones(1,6)*x(1:6)-(r(2)-p(2))*ones(1,2)*x(7:8)-(r(3)-p(3))*x(9);
-    prob.Constraints.con1=p(1)*ones(1,6)*x(1:6)+p(2)*ones(1,2)*x(7:8)+p(3)*x(9)+dp<=total(i);%总费用
+    prob.Constraints.con1=p(1)*ones(1,6)*x(1:6)+p(2)*ones(1,2)*x(7:8)+p(3)*x(9)+dp<=total(i);%总费用
-    prob.Constraints.con2=ot(1,1)*ones(1,3)*x(1:3)+ot(1,2)*x(7)<=6000;%A1
+    prob.Constraints.con2=ot(1,1)*ones(1,3)*x(1:3)+ot(1,2)*x(7)<=6000;%A1
-    prob.Constraints.con3=ot(2,1)*ones(1,3)*x(4:6)+ot(2,2)*x(8)+ot(2,3)*x(9)<=10000;
+    prob.Constraints.con3=ot(2,1)*ones(1,3)*x(4:6)+ot(2,2)*x(8)+ot(2,3)*x(9)<=10000;
-    prob.Constraints.con4=ot(3,1)*(x(1)+x(4))+ot(3,2)*(x(7)+x(8))<=4000;
+    prob.Constraints.con4=ot(3,1)*(x(1)+x(4))+ot(3,2)*(x(7)+x(8))<=4000;
-    prob.Constraints.con5=ot(4,1)*(x(2)+x(5))+ot(4,3)*x(9)<=7000;
+    prob.Constraints.con5=ot(4,1)*(x(2)+x(5))+ot(4,3)*x(9)<=7000;
-    prob.Constraints.con6=ot(5,1)*(x(3)+x(6))<=4000;
+    prob.Constraints.con6=ot(5,1)*(x(3)+x(6))<=4000;
-    [sol,fval,flag,out]=solve(prob);
+    [sol,fval,flag,out]=solve(prob);
-    xx=sol.x;
+    xx=sol.x;
-
+
-    X=[X,xx];
+    X=[X,xx];
-    Q=[Q,-fval];
+    Q=[Q,-fval];
-
+
-end
+end
-plot(total,Q,'*-');
+plot(total,Q,'*-');
-
+
-
+
-
+
-
+
-
+

姚安欣/Mengtkl.m

@@ -1,61 +1,61 @@
-
+
-%27页例题
+%27页例题
-clc; clear; close all;
+clc; clear; close all;
-n = 10000; % 使用较小的 n 值以便更容易可视化
+n = 10000; % 使用较小的 n 值以便更容易可视化
-x = unifrnd(0, 12, [1, n]);
+x = unifrnd(0, 12, [1, n]);
-y = unifrnd(0, 9, [1, n]);
+y = unifrnd(0, 9, [1, n]);
-ans=sum(y < x.^2 & x <= 3)+sum(y < 12 - x & x >= 3);
+ans=sum(y < x.^2 & x <= 3)+sum(y < 12 - x & x >= 3);
-ans=ans/n;
+ans=ans/n;
-% 找出满足条件的点
+% 找出满足条件的点
-condition1 = y < x.^2 & x <= 3;
+condition1 = y < x.^2 & x <= 3;
-condition2 = y < 12 - x & x >= 3;
+condition2 = y < 12 - x & x >= 3;
-condition_met = condition1 | condition2; % 满足任一条件的点
+condition_met = condition1 | condition2; % 满足任一条件的点
-condition_not_met = ~condition_met; % 不满足任何条件的点
+condition_not_met = ~condition_met; % 不满足任何条件的点
-
+
-% 创建图形窗口
+% 创建图形窗口
-figure;
+figure;
-hold on;%在同一张图上绘图
+hold on;%在同一张图上绘图
-
+
-% 绘制不满足任何条件的点
+% 绘制不满足任何条件的点
-scatter(x(condition_not_met), y(condition_not_met), 'k.'); %  k----黑色  .----绘制样式
+scatter(x(condition_not_met), y(condition_not_met), 'k.'); %  k----黑色  .----绘制样式
-%scatter绘制散点图
+%scatter绘制散点图
-%x(condition_not_met) 会返回一个新的向量，其中只包含 x 中对应 condition_not_met 为 true 的元素。
+%x(condition_not_met) 会返回一个新的向量，其中只包含 x 中对应 condition_not_met 为 true 的元素。
-
+
-% 绘制满足第一个条件的点
+% 绘制满足第一个条件的点
-scatter(x(condition1), y(condition1), 'r.'); % 红色
+scatter(x(condition1), y(condition1), 'r.'); % 红色
-
+
-% 绘制满足第二个条件的点
+% 绘制满足第二个条件的点
-scatter(x(condition2), y(condition2), 'b.'); % 蓝色
+scatter(x(condition2), y(condition2), 'b.'); % 蓝色
-
+
-% 添加图例和标签
+% 添加图例和标签
-legend('不满足任何条件的点', '满足 y < x^2 且 x <= 3 的点', '满足 y < 12 - x 且 x >= 3 的点');
+legend('不满足任何条件的点', '满足 y < x^2 且 x <= 3 的点', '满足 y < 12 - x 且 x >= 3 的点');
-xlabel('x');
+xlabel('x');
-ylabel('y');
+ylabel('y');
-title('随机生成的点和满足条件的点');
+title('随机生成的点和满足条件的点');
-hold off;
+hold off;
-
+
-
+
-%蒙特卡洛法求圆周率qw
+%蒙特卡洛法求圆周率qw
-clc;clear;close all;
+clc;clear;close all;
-
+
-n=10^5;
+n=10^5;
-x=unifrnd(-1,1,[1,n]);
+x=unifrnd(-1,1,[1,n]);
-y=unifrnd(-1,1,[1,n]);
+y=unifrnd(-1,1,[1,n]);
-con1=x.^2+y.^2<=1;
+con1=x.^2+y.^2<=1;
-con2=~con1;
+con2=~con1;
-ans=sum(x.^2+y.^2<=1);
+ans=sum(x.^2+y.^2<=1);
-ans=ans/n*4;
+ans=ans/n*4;
-
+
-figure ;
+figure ;
-hold on;
+hold on;
-
+
-scatter(x(con1),y(con1),'r.');
+scatter(x(con1),y(con1),'r.');
-scatter(x(con2),y(con2),'k.');
+scatter(x(con2),y(con2),'k.');
-
+
-hold off;
+hold off;
-
+
-
+
-
+
-
+
-
+
-
+

姚安欣/NLP.m

@@ -1,135 +1,135 @@
-
+
-%非线性规划
+%非线性规划
-
+
-%P42例题
+%P42例题
-
+
-%求解
+%求解
-clc,clear,format long g,close all;
+clc,clear,format long g,close all;
-syms x1 x2;%符号变量(用于定义函数、求解方程、计算导数和积分)
+syms x1 x2;%符号变量(用于定义函数、求解方程、计算导数和积分)
-f=(339-0.01*x1-0.003*x2)*x1+(399-0.004*x1-0.01*x2)*x2-(195*x1+225*x2+400000);
+f=(339-0.01*x1-0.003*x2)*x1+(399-0.004*x1-0.01*x2)*x2-(195*x1+225*x2+400000);
-f=simplify(f);
+f=simplify(f);
-f1=diff(f,x1),f2=diff(f,x2);
+f1=diff(f,x1),f2=diff(f,x2);
-[x10,x20]=solve(f1,f2);%求驻点即使f1,f2为0的点
+[x10,x20]=solve(f1,f2);%求驻点即使f1,f2为0的点
-x10=round(double(x10));%四舍五入
+x10=round(double(x10));%四舍五入
-x20=round(double(x20));
+x20=round(double(x20));
-f0=subs(f,{x1,x2},{x10,x20});
+f0=subs(f,{x1,x2},{x10,x20});
-f0=double(f0);
+f0=double(f0);
-subplot(121);
+subplot(121);
-fmesh(f,[0,10000,0,10000]);%3维
+fmesh(f,[0,10000,0,10000]);%3维
-%xlabel('x1','Interpreter','latex');
+%xlabel('x1','Interpreter','latex');
-%ylabel('y1','Interpreter','latex');
+%ylabel('y1','Interpreter','latex');
-xlabel('x1');
+xlabel('x1');
-ylabel('y1');
+ylabel('y1');
-subplot(122);
+subplot(122);
-L=fcontour(f,[0,10000,0,10000]);
+L=fcontour(f,[0,10000,0,10000]);
-contour(L.XData,L.YData,L.ZData,'ShowText','on');
+contour(L.XData,L.YData,L.ZData,'ShowText','on');
-xlabel('x1');
+xlabel('x1');
-ylabel('y1');
+ylabel('y1');
-p1=339-0.01*x10-0.003*x20;
+p1=339-0.01*x10-0.003*x20;
-p2=399-0.004*x10-0.01*x20;
+p2=399-0.004*x10-0.01*x20;
-c=195*x10+225*x20+400000;
+c=195*x10+225*x20+400000;
-rate=f0/c;
+rate=f0/c;
-
+
-%灵敏性检验 分析 最优解 对 19村彩电价格下降幅度 的灵敏性
+%灵敏性检验 分析 最优解 对 19村彩电价格下降幅度 的灵敏性
-clc,clear,close all,format long g;
+clc,clear,close all,format long g;
-syms x1 x2 a;
+syms x1 x2 a;
-X1=4735,X2=7043,AA=0.01,FF=553641;
+X1=4735,X2=7043,AA=0.01,FF=553641;
-f=(339-a*x1-0.003*x2)*x1+(399-0.004*x1-0.01*x2)*x2-(195*x1+225*x2+400000);
+f=(339-a*x1-0.003*x2)*x1+(399-0.004*x1-0.01*x2)*x2-(195*x1+225*x2+400000);
-f1=diff(f,x1),f2=diff(f,x2);
+f1=diff(f,x1),f2=diff(f,x2);
-[x10,x20]=solve(f1,f2);
+[x10,x20]=solve(f1,f2);
-
+
-subplot(121);
+subplot(121);
-fplot(x10,[0.002,0.02]);
+fplot(x10,[0.002,0.02]);
-xlabel('a');
+xlabel('a');
-ylabel('x1');
+ylabel('x1');
-
+
-subplot(122);
+subplot(122);
-fplot(x20,[0.002,0.02]);
+fplot(x20,[0.002,0.02]);
-xlabel('a');
+xlabel('a');
-ylabel('x2');
+ylabel('x2');
-
+
-dx1=diff(x10,a);
+dx1=diff(x10,a);
-dx10=subs(dx1,a,0.01),dx10=double(dx10);
+dx10=subs(dx1,a,0.01),dx10=double(dx10);
-sx1a=dx10*AA/X1;
+sx1a=dx10*AA/X1;
-
+
-dx2=diff(x20,a);
+dx2=diff(x20,a);
-dx20=subs(dx2,a,0.01),dx20=double(dx20);
+dx20=subs(dx2,a,0.01),dx20=double(dx20);
-sx2a=dx20*AA/X2;
+sx2a=dx20*AA/X2;
-
+
-F=subs(f,[x1,x2],[x10,x20]);
+F=subs(f,[x1,x2],[x10,x20]);
-F=simplify(F);
+F=simplify(F);
-figure,fplot(F,[0.002,0.02]);
+figure,fplot(F,[0.002,0.02]);
-xlabel('a');
+xlabel('a');
-ylabel('profit');
+ylabel('profit');
-dpro=diff(F,a);
+dpro=diff(F,a);
-dpro0=subs(F,a,0.01);
+dpro0=subs(F,a,0.01);
-sya=dpro0*AA/FF;
+sya=dpro0*AA/FF;
-
+
-%当a提高10%
+%当a提高10%
-f3=subs(f,[x1,x2,a],[X1,X2,0.011]),f3=double(f3);
+f3=subs(f,[x1,x2,a],[X1,X2,0.011]),f3=double(f3);
-f4=subs(F,[a],0.011);
+f4=subs(F,[a],0.011);
-%相对误差
+%相对误差
-delta=(f4-f3)/f4;
+delta=(f4-f3)/f4;
-delta=double(delta);
+delta=double(delta);
-%相对误差很小，说明在19寸的价格下降幅度在一定范围内发生变化时，仍然按照不变时的最优解来确定生产方案，仍然可以获得几乎最优的利润。
+%相对误差很小，说明在19寸的价格下降幅度在一定范围内发生变化时，仍然按照不变时的最优解来确定生产方案，仍然可以获得几乎最优的利润。
-
+
-
+
-
+
-
+
-
+
-%使用fmincon的简单例子
+%使用fmincon的简单例子
-
+
-%蒙特卡罗确定初始值
+%蒙特卡罗确定初始值
-clc,clear,close all;
+clc,clear,close all;
-n=10^7;
+n=10^7;
-x1=unifrnd(0,12,[1,n]);
+x1=unifrnd(0,12,[1,n]);
-x2=unifrnd(0,12,[1,n]);
+x2=unifrnd(0,12,[1,n]);
-fmin=+inf;
+fmin=+inf;
-for i=1:n
+for i=1:n
-    x=[x1(i),x2(i)];
+    x=[x1(i),x2(i)];
-    if ((x1(i)-1)^2-x2(i)-6<=0) & (-2*x1(i)+3*x2(i)-6<=0)
+    if ((x1(i)-1)^2-x2(i)-6<=0) & (-2*x1(i)+3*x2(i)-6<=0)
-        temp=x(1)^2+x(2)^2-x(1)*x(2)-2*x(1)-5*x(2);
+        temp=x(1)^2+x(2)^2-x(1)*x(2)-2*x(1)-5*x(2);
-        if temp<fmin
+        if temp<fmin
-           fmin=temp;
+           fmin=temp;
-           x0=x;
+           x0=x;
-        end
+        end
-    end
+    end
-end
+end
-
+
-%[x,fval]= fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options)
+%[x,fval]= fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options)
-
+
-A=[-2,3];
+A=[-2,3];
-b=6;
+b=6;
-[x,fval]=fmincon(@fun1,x0,A,b,[],[],[],[],@nonlfun1);
+[x,fval]=fmincon(@fun1,x0,A,b,[],[],[],[],@nonlfun1);
-
+
-%option=optimoptions('fmincon','Algorithm','sqp');
+%option=optimoptions('fmincon','Algorithm','sqp');
-%[x,fval]=fmincon(@fun1,x0,A,b,[],[],[],[],@nonlfun1,option);
+%[x,fval]=fmincon(@fun1,x0,A,b,[],[],[],[],@nonlfun1,option);
-
+
-disp('Thw optimal solution is:');
+disp('Thw optimal solution is:');
-disp(x);                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                        
+disp(x);                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                        
-disp('The minimum value is:');
+disp('The minimum value is:');
-disp(fval);
+disp(fval);
-figure;
+figure;
-fmesh(@(x1, x2) arrayfun(@(x1, x2) fun1([x1, x2]), x1, x2), [-4, 4, -4, 4]);%用法
+fmesh(@(x1, x2) arrayfun(@(x1, x2) fun1([x1, x2]), x1, x2), [-4, 4, -4, 4]);%用法
-
+
-xlabel('x1');
+xlabel('x1');
-ylabel('x2');
+ylabel('x2');
-zlabel('f(x1,x2)');
+zlabel('f(x1,x2)');
-
+
-
+
-%非线性规划---料场选址
+%非线性规划---料场选址
-clc,clear,close all;
+clc,clear,close all;
-a=[1.25,8.75,0.5,5.75,3,7.25];
+a=[1.25,8.75,0.5,5.75,3,7.25];
-b=[1.25,0.75,4.75,5,6.5,7.25];
+b=[1.25,0.75,4.75,5,6.5,7.25];
-A=[1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0;0 0 0 0 0 0 1 1 1 1 1 1 0 0 0 0];
+A=[1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0;0 0 0 0 0 0 1 1 1 1 1 1 0 0 0 0];
-bb=[20;20];
+bb=[20;20];
-Aeq=[eye(6),eye(6),zeros(6,4)];
+Aeq=[eye(6),eye(6),zeros(6,4)];
-beq=[3;5;4;7;6;11;];
+beq=[3;5;4;7;6;11;];
-lb=[zeros(12,1);-inf;-inf;-inf;-inf];
+lb=[zeros(12,1);-inf;-inf;-inf;-inf];
-x0=[3,5,0,7,0,1,0,0,4,0,6,10,5,1,2,7];
+x0=[3,5,0,7,0,1,0,0,4,0,6,10,5,1,2,7];
-[xans,fvans]=fmincon(@fun2,x0,A,bb,Aeq,beq,lb,[]);
+[xans,fvans]=fmincon(@fun2,x0,A,bb,Aeq,beq,lb,[]);
-
+
-
+
-
+
-
+
-
+

姚安欣/PCA.m

@@ -0,0 +1,96 @@
+clc,clear,format short,close all;
+load data_mh.mat;
+[n,p]=size(x);
+
+%标准化
+X=zscore(x);
+
+%协方差矩阵/相关系数矩阵
+R=cov(X);
+[V,D]=eig(R);%[特征向量，特征值]
+
+lambda=diag(D);
+lambda=lambda(end:-1:1);
+
+total_contri=sum(lambda);
+
+cum_contri=cumsum(lambda);
+
+contri_rate=cum_contri/total_contri;
+
+V1=rot90(V)';
+disp(V1);
+
+c1=V1(:,1);
+c2=V1(:,2);
+
+XX=X*c1;
+YY=X*c2;
+
+
+figure;
+scatter(XX, YY, 'filled');
+xlabel('第一主成分');
+ylabel('第二主成分');
+title('主成分得分图');
+grid on;
+
+
+
+
+
+
+%散点
+figure;
+scatter(XX,YY,'k.');
+
+
+figure;
+plot(c1, c2, 'o', 'MarkerSize', 8, 'MarkerFaceColor', 'b');
+text(c1, c2, cellstr(num2str((1:p)')), 'VerticalAlignment','bottom', 'HorizontalAlignment','right');
+xlabel('第一主成分');
+ylabel('第二主成分');
+title('主成分载荷图');
+grid on;
+
+
+
+
+
+
+
+
+
+
+
+
+%主成分得分图（Score Plot）
+X = randn(100, 5); 
+
+% 主成分分析
+[coeff, score, latent, tsquared, explained, mu] = pca(X);
+
+% 绘制前两个主成分的得分图
+figure;
+scatter(score(:,1), score(:,2));
+xlabel('第一主成分');
+ylabel('第二主成分');
+title('主成分得分图');
+grid on;
+
+% 添加数据点标签
+text(score(:,1), score(:,2), num2str((1:size(score,1))'), 'FontSize', 8);
+
+
+
+
+
+
+
+
+
+
+
+
+
+

姚安欣/factor.m

@@ -0,0 +1,46 @@
+clc,clear,close all;
+
+%相关系数矩阵
+r=[ 1.000,0.577,0.509,0.387,0.462
+    0.577,1.000,1.599,0.389,0.322
+    0.509,0.599,1.000,0.436,0.426
+    0.387,0.389,0.436,1.000,0.523
+    0.462,0.322,0.426,0.523,1.000];
+
+[vec1,val,rate]=pcacov(r);%特征向量、特征值、贡献率
+f1=repmat(sign(sum(vec1)),size(vec1,1),1);%调整符号
+vec2=vec1.*f1;%是用  .*
+f2=repmat(sqrt(val)',size(vec2,1),1);
+a=vec2.*f2;%载荷矩阵
+a1=a(:,1);
+tcha1=diag(r-a1*a1');
+a2=a(:,[1,2]);
+tcha2=diag(r-a2*a2');
+ccha2=r-a2*a2'-diag(tcha2);
+con=cumsum(rate);
+
+
+
+
+clc,clear,close all;
+load data_mh.mat;
+[n,p]=size(x);
+%标准化
+X=zscore(x);
+
+%相关系数矩阵
+r=cov(X);
+
+[vec1,val,rate]=pcacov(r);%特征向量、特征值、贡献率
+f1=repmat(sign(sum(vec1)),size(vec1,1),1);%调整符号
+vec2=vec1.*f1;%是用  .*
+f2=repmat(sqrt(val)',size(vec2,1),1);
+a=vec2.*f2;%载荷矩阵
+a1=a(:,1);
+tcha1=diag(r-a1*a1');
+a2=a(:,[1,2]);
+tcha2=diag(r-a2*a2');
+ccha2=r-a2*a2'-diag(tcha2);
+con=cumsum(rate);
+
+

姚安欣/fun1.m

@@ -1,3 +1,3 @@
-function f=fun1(x)
+function f=fun1(x)
-    f=x(1)*x(1)+x(2)*x(2)-x(1)*x(2)-2*x(1)-5*x(2);
+    f=x(1)*x(1)+x(2)*x(2)-x(1)*x(2)-2*x(1)-5*x(2);
-end
+end

姚安欣/fun2.m

@@ -1,15 +1,15 @@
-function f=fun2(x)
+function f=fun2(x)
-a=[1.25,8.75,0.5,5.75,3,7.25];
+a=[1.25,8.75,0.5,5.75,3,7.25];
-b=[1.25,0.75,4.75,5,6.5,7.25];
+b=[1.25,0.75,4.75,5,6.5,7.25];
- f1 = 0;  % 初始化f1  
+ f1 = 0;  % 初始化f1  
- f2 = 0;  % 初始化f2  
+ f2 = 0;  % 初始化f2  
-for i=1:6
+for i=1:6
-    s=sqrt((a(i)-x(13))^2+(b(i)-x(14))^2);
+    s=sqrt((a(i)-x(13))^2+(b(i)-x(14))^2);
-    f1=x(i)*s+f1;
+    f1=x(i)*s+f1;
-end
+end
-for i=7:12
+for i=7:12
-    s=sqrt((a(i-6)-x(15))^2+(b(i-6)-x(16))^2);
+    s=sqrt((a(i-6)-x(15))^2+(b(i-6)-x(16))^2);
-    f2=s*x(i)+f2;
+    f2=s*x(i)+f2;
-end
+end
-f=f1+f2;
+f=f1+f2;
 end

姚安欣/intLP.m

@@ -1,49 +1,49 @@
-%整数规划
+%整数规划
-%P23 ej2.5
+%P23 ej2.5
-%optimproblem解法
+%optimproblem解法
-clc; clear; close all; format long g;
+clc; clear; close all; format long g;
-prob=optimproblem;
+prob=optimproblem;
-x=optimvar('x',6,'LowerBound',0,'Type','integer');
+x=optimvar('x',6,'LowerBound',0,'Type','integer');
-prob.Objective=sum(x);
+prob.Objective=sum(x);
-cnt=[35,40,50,45,55,30];
+cnt=[35,40,50,45,55,30];
-con=optimconstr(6);
+con=optimconstr(6);
-con(1)=x(1)+x(6)>=35;
+con(1)=x(1)+x(6)>=35;
-for i=1:5
+for i=1:5
-    con(i+1)=x(i)+x(i+1)>=cnt(i+1);
+    con(i+1)=x(i)+x(i+1)>=cnt(i+1);
-end
+end
-prob.Constraints.con=con;
+prob.Constraints.con=con;
-[sol,fval,flag,out]=solve(prob);
+[sol,fval,flag,out]=solve(prob);
-X=sol.x;
+X=sol.x;
-
+
-
+
-%linprog解法
+%linprog解法
-clc; clear; close all; format long g;
+clc; clear; close all; format long g;
-f=[1,1,1,1,1,1];
+f=[1,1,1,1,1,1];
-intcon=[1:6];
+intcon=[1:6];
-A=zeros(6,6);
+A=zeros(6,6);
-A(1,1)=-1;
+A(1,1)=-1;
-A(1,6)=-1;
+A(1,6)=-1;
-for i=1:5
+for i=1:5
-    A(i+1,i)=-1;
+    A(i+1,i)=-1;
-    A(i+1,i+1)=-1;
+    A(i+1,i+1)=-1;
-end
+end
-lb=zeros(6,1);%注意不是lb=0
+lb=zeros(6,1);%注意不是lb=0
-b=[-35;-40;-50;-45;-55;-30];
+b=[-35;-40;-50;-45;-55;-30];
-[x,fval]=intlinprog(f,intcon,A,b,[],[],lb);
+[x,fval]=intlinprog(f,intcon,A,b,[],[],lb);
-
+
-
+
-
+
-%01规划
+%01规划
-%背包问题
+%背包问题
-
+
-f=-[540,200,180,350,60,150,280,450,320,120];
+f=-[540,200,180,350,60,150,280,450,320,120];
-intcon=1:10;
+intcon=1:10;
-lb=zeros(10);
+lb=zeros(10);
-ub=ones(10);
+ub=ones(10);
-A=[6,3,4,5,1,2,3,5,4,2];
+A=[6,3,4,5,1,2,3,5,4,2];
-b=30;
+b=30;
-[x,fval]=intlinprog(f,intcon,A,b,[],[],lb,ub);
+[x,fval]=intlinprog(f,intcon,A,b,[],[],lb,ub);
-
+
-
+
-
+
-
+

姚安欣/mm1.m

@@ -1,140 +1,140 @@
-%投资利益与风险1998A题
+%投资利益与风险1998A题
-
+
-
+
-%模型一 给定风险承受程度，求最大利益
+%模型一 给定风险承受程度，求最大利益
-f=[-0.05,-0.27,-0.19,-0.185,-0.185];
+f=[-0.05,-0.27,-0.19,-0.185,-0.185];
-
+
-%A矩阵
+%A矩阵
-%A=[0,0.025,0.015,0.055,0.026];         %错误
+%A=[0,0.025,0.015,0.055,0.026];         %错误
-%b=[1,1,1,1,1];
+%b=[1,1,1,1,1];
-A=[zeros(4,1),diag([0.025,0.015,0.055,0.026])];%不等式约束条件矩阵
+A=[zeros(4,1),diag([0.025,0.015,0.055,0.026])];%不等式约束条件矩阵
-
+
-%Aeq、beq
+%Aeq、beq
-Aeq=[1,1.01,1.02,1.045,1.065];
+Aeq=[1,1.01,1.02,1.045,1.065];
-beq=1;
+beq=1;
-
+
-%lb
+%lb
-%lb=0;                  %错误
+%lb=0;                  %错误
-lb=zeros(5,1);
+lb=zeros(5,1);
-
+
-
+
-%可承担风险率a
+%可承担风险率a
-a=(0:0.001:0.05);
+a=(0:0.001:0.05);
-
+
-%保存最优解
+%保存最优解
-Q=zeros(1,length(a));
+Q=zeros(1,length(a));
-xx=[];%空矩阵存放最优解对应x的值
+xx=[];%空矩阵存放最优解对应x的值
-for i=1:length(a)
+for i=1:length(a)
-    b=a(i)*ones(4,1);
+    b=a(i)*ones(4,1);
-    [x,y]=linprog(f,A,b,Aeq,beq,lb);
+    [x,y]=linprog(f,A,b,Aeq,beq,lb);
-    Q(i)=-y;%注意取负！！！
+    Q(i)=-y;%注意取负！！！
-    xx=[xx;x'];
+    xx=[xx;x'];
-end
+end
-plot(a,Q,'*r');
+plot(a,Q,'*r');
-xlabel("风险率");
+xlabel("风险率");
-ylabel("最大收益");
+ylabel("最大收益");
-
+
-
+
-%模型二 收益、风险按权重组合
+%模型二 收益、风险按权重组合
-% f0=[-0.05,-0.27,-0.19,-0.185,-0.185];
+% f0=[-0.05,-0.27,-0.19,-0.185,-0.185];
-% w=(0:0.1:1);
+% w=(0:0.1:1);
-% Aeq=[1,1.01,1.02,1.045,1.065,0];
+% Aeq=[1,1.01,1.02,1.045,1.065,0];
-% beq=1;
+% beq=1;
-% lb=0;
+% lb=0;
-% xx=[];
+% xx=[];
-% Q=zeros(1,length(w));
+% Q=zeros(1,length(w));
-% A=[zeros(5,1),diag([0.025,0.025,0.055,0.065,0])];
+% A=[zeros(5,1),diag([0.025,0.025,0.055,0.065,0])];
-% b=ones(5,1);
+% b=ones(5,1);
-% for i=1:length(w)
+% for i=1:length(w)
-%     f=[-w(i)*f0,1-w(i)];
+%     f=[-w(i)*f0,1-w(i)];
-%     b=x(end)*b;
+%     b=x(end)*b;
-%     [x,y]=linprog(f,A,b,Aeq,beq,lb);
+%     [x,y]=linprog(f,A,b,Aeq,beq,lb);
-%     Q(i)=-y;
+%     Q(i)=-y;
-%     xx=[xx,x'];
+%     xx=[xx,x'];
-% end
+% end
-% plot(w,Q,'*r');
+% plot(w,Q,'*r');
-
+
-%模型二 收益、风险按权重组合
+%模型二 收益、风险按权重组合
-clc; clear; close all; format long g;
+clc; clear; close all; format long g;
-
+
-M = 10000;
+M = 10000;
-prob = optimproblem;
+prob = optimproblem;
-x = optimvar('x', 6, 1, 'LowerBound', 0);
+x = optimvar('x', 6, 1, 'LowerBound', 0);
-r = [0.05, 0.28, 0.21, 0.23, 0.25];
+r = [0.05, 0.28, 0.21, 0.23, 0.25];
-p = [0, 0.01, 0.02, 0.045, 0.065];
+p = [0, 0.01, 0.02, 0.045, 0.065];
-q = [0, 0.025, 0.015, 0.055, 0.026]';
+q = [0, 0.025, 0.015, 0.055, 0.026]';
-%w = 0:0.1:1;
+%w = 0:0.1:1;
-w = 0.7:0.03:1;
+w = 0.7:0.03:1;
-V = [];
+V = [];
-Q = [];
+Q = [];
-X = [];
+X = [];
-
+
-prob.Constraints.con1 = (1 + p) * x(1:end-1) == M;
+prob.Constraints.con1 = (1 + p) * x(1:end-1) == M;
-prob.Constraints.con2 = (q(2:end).* x(2:end-1))<= x(end);        %下标从1开始
+prob.Constraints.con2 = (q(2:end).* x(2:end-1))<= x(end);        %下标从1开始
-
+
-for i = 1:length(w)
+for i = 1:length(w)
-    prob.Objective = w(i) * x(end) - (1 - w(i)) * (r - p) * x(1:end-1);     %注意大小写
+    prob.Objective = w(i) * x(end) - (1 - w(i)) * (r - p) * x(1:end-1);     %注意大小写
-    [sol, fval, flag, out] = solve(prob);
+    [sol, fval, flag, out] = solve(prob);
-    xx = sol.x;
+    xx = sol.x;
-    V = [V, max(q.* xx(1:end-1))];
+    V = [V, max(q.* xx(1:end-1))];
-    Q = [Q, (r - p) * xx(1:end-1)];
+    Q = [Q, (r - p) * xx(1:end-1)];
-    X = [X, xx];
+    X = [X, xx];
-    
+    
-    plot(V, Q, '*-');
+    plot(V, Q, '*-');
-    grid on;
+    grid on;
-    xlabel('风险');
+    xlabel('风险');
-    ylabel('收益');
+    ylabel('收益');
-end
+end
-
+
-
+
-
+
-%模型三：达到一定盈利水平，极小化风险
+%模型三：达到一定盈利水平，极小化风险
-
+
-clc; clear; close all; format long g;
+clc; clear; close all; format long g;
-M=10000;
+M=10000;
-k=1500:100:3000;
+k=1500:100:3000;
-prob = optimproblem;
+prob = optimproblem;
-x = optimvar('x', 5, 1, 'LowerBound', 0);%下界为0
+x = optimvar('x', 5, 1, 'LowerBound', 0);%下界为0
-r = [0.05, 0.28, 0.21, 0.23, 0.25];
+r = [0.05, 0.28, 0.21, 0.23, 0.25];
-p = [0, 0.01, 0.02, 0.045, 0.065];
+p = [0, 0.01, 0.02, 0.045, 0.065];
-q = [0, 0.025, 0.015, 0.055, 0.026]';
+q = [0, 0.025, 0.015, 0.055, 0.026]';
-
+
-V = [];
+V = [];
-Q = [];
+Q = [];
-X = [];
+X = [];
-t = optimvar('t', 'LowerBound', 0);
+t = optimvar('t', 'LowerBound', 0);
-
+
-%prob.Objective=max(q.*x);%极小化风险
+%prob.Objective=max(q.*x);%极小化风险
-for i=1:length(k)
+for i=1:length(k)
-    prob.Objective = t; % 极小化风险
+    prob.Objective = t; % 极小化风险
-    prob.Constraints.con1 = (1 + p) * x == M;
+    prob.Constraints.con1 = (1 + p) * x == M;
-    prob.Constraints.con2=((r-p)*x>=k(i));%达到一定盈利水平
+    prob.Constraints.con2=((r-p)*x>=k(i));%达到一定盈利水平
-    prob.Constraints.con3=(q.*x<=t);
+    prob.Constraints.con3=(q.*x<=t);
-    [sol,fval,flag,out]=solve(prob);
+    [sol,fval,flag,out]=solve(prob);
-    if flag==1
+    if flag==1
-        xx=sol.x;
+        xx=sol.x;
-        X=[X,xx];
+        X=[X,xx];
-        Q=[Q,(r-p)*xx];
+        Q=[Q,(r-p)*xx];
-        V=[V,fval];
+        V=[V,fval];
-    else
+    else
-        xx=-1*ones(5,1);
+        xx=-1*ones(5,1);
-        X=[X,xx];
+        X=[X,xx];
-        Q=[Q,-1];
+        Q=[Q,-1];
-        V=[V,-1];
+        V=[V,-1];
-    end
+    end
-
+
-
+
-end
+end
-plot(k,V,'*-');
+plot(k,V,'*-');
-
+
-
+
-
+
-
+
-
+
-
+
-
+
-
+
-
+
-
+
-
+
-
+

姚安欣/mm2_graph.m

@@ -1,94 +1,94 @@
-%ej.4.10仓库选址
+%ej.4.10仓库选址
-clc; clear; close all; format long g;
+clc; clear; close all; format long g;
-a=zeros(6);
+a=zeros(6);
-a(1,[2,5])=[20,15];
+a(1,[2,5])=[20,15];
-a(2,[1,3,4,5])=[20,20,60,25];
+a(2,[1,3,4,5])=[20,20,60,25];
-a(3,[2,4,5])=[20,30,18];
+a(3,[2,4,5])=[20,30,18];
-a(4,[2,3])=[60,30];
+a(4,[2,3])=[60,30];
-a(5,[1,2,3,6])=[15,25,18,15];
+a(5,[1,2,3,6])=[15,25,18,15];
-a(6,5)=15;
+a(6,5)=15;
-s = cellstr(strcat('销售点',int2str([1:6]')));
+s = cellstr(strcat('销售点',int2str([1:6]')));
-G=graph(a,s);
+G=graph(a,s);
-P = plot(G,'layout','force','EdgeColor','k','NodeFontSize',12);
+P = plot(G,'layout','force','EdgeColor','k','NodeFontSize',12);
-D=[];
+D=[];
-for i=1:6
+for i=1:6
-    temp=0;
+    temp=0;
-    for j=1:6
+    for j=1:6
-        [path,d]=shortestpath(G,i,j);
+        [path,d]=shortestpath(G,i,j);
-        if d>temp;
+        if d>temp;
-            temp=d;
+            temp=d;
-        end
+        end
-    end
+    end
-    D=[D,temp];
+    D=[D,temp];
-end
+end
-[minVal,minInd]=min(D);
+[minVal,minInd]=min(D);
-disp('仓库应建在销售点');
+disp('仓库应建在销售点');
-disp(minInd);
+disp(minInd);
-
+
-
+
-%答案版
+%答案版
-clc; clear; close all;
+clc; clear; close all;
-
+
-% 初始化邻接矩阵
+% 初始化邻接矩阵
-a = zeros(6);
+a = zeros(6);
-a(1, [2, 5]) = [20, 15];
+a(1, [2, 5]) = [20, 15];
-a(2, 3:5) = [20, 60, 25];
+a(2, 3:5) = [20, 60, 25];
-a(3, [4, 5]) = [30, 18];
+a(3, [4, 5]) = [30, 18];
-a(5, 6) = 15;
+a(5, 6) = 15;
-
+
-% 创建节点标签
+% 创建节点标签
-s = cellstr(strcat('V', int2str([1:6]')));
+s = cellstr(strcat('V', int2str([1:6]')));
-
+
-% 创建图
+% 创建图
-G = graph(a, s, 'upper');
+G = graph(a, s, 'upper');
-
+
-% 计算所有节点之间的最短路径距离矩阵
+% 计算所有节点之间的最短路径距离矩阵
-d = distances(G);
+d = distances(G);
-
+
-% 计算从每个节点到其他节点的最大最短路径距离
+% 计算从每个节点到其他节点的最大最短路径距离
-D = max(d, [], 2)';     %第三个参数 2 指定沿着矩阵的第二个维度（即列）进行操作。
+D = max(d, [], 2)';     %第三个参数 2 指定沿着矩阵的第二个维度（即列）进行操作。
-
+
-% 找到最大最短路径距离中的最小值及其下标
+% 找到最大最短路径距离中的最小值及其下标
-[minD, minIndex] = min(D);
+[minD, minIndex] = min(D);
-
+
-% 绘制图
+% 绘制图
-plot(G, 'EdgeLabel', G.Edges.Weight, 'Layout', 'force');%设置边的标签
+plot(G, 'EdgeLabel', G.Edges.Weight, 'Layout', 'force');%设置边的标签
-
+
-% 显示结果
+% 显示结果
-fprintf('从每个节点到任意其他节点的最大最短路径距离: \n');
+fprintf('从每个节点到任意其他节点的最大最短路径距离: \n');
-disp(D);
+disp(D);
-fprintf('最大最短路径距离中的最小值为: %f\n', minD);
+fprintf('最大最短路径距离中的最小值为: %f\n', minD);
-fprintf('对应的节点为: V%d\n', minIndex);
+fprintf('对应的节点为: V%d\n', minIndex);
-
+
-%设备更新问题；
+%设备更新问题；
-clc;clear;close all;
+clc;clear;close all;
-a=zeros(6);
+a=zeros(6);
-a(1,[2:6])=[15,20,27,37,54];
+a(1,[2:6])=[15,20,27,37,54];
-a(2,[3:6])=[15,20,27,37];
+a(2,[3:6])=[15,20,27,37];
-a(3,[4:6])=[16,21,28];
+a(3,[4:6])=[16,21,28];
-a(4,[5:6])=[16,21];
+a(4,[5:6])=[16,21];
-a(5,6)=17;
+a(5,6)=17;
-s=cellstr(strcat(int2str((1:6)'),'年初'));
+s=cellstr(strcat(int2str((1:6)'),'年初'));
-G=digraph(a,s);
+G=digraph(a,s);
-P = plot(G,'layout','force','EdgeColor','k','NodeFontSize',12);
+P = plot(G,'layout','force','EdgeColor','k','NodeFontSize',12);
-[path,d]=shortestpath(G,1,6);%注意参数个数
+[path,d]=shortestpath(G,1,6);%注意参数个数
-highlight(P,path,"EdgeColor","red",'LineWidth',3.5);
+highlight(P,path,"EdgeColor","red",'LineWidth',3.5);
-
+
-
+
-%最小生成树
+%最小生成树
-clc;clear;
+clc;clear;
-a=zeros(9);
+a=zeros(9);
-a(1,[2:9])=[2 1 3 4 4 2 5 4];       % 顶点1到其他顶点的边的权重
+a(1,[2:9])=[2 1 3 4 4 2 5 4];       % 顶点1到其他顶点的边的权重
-a(2,[3 9])=[4 1];        % 顶点2到顶点3、顶点9的边的权重
+a(2,[3 9])=[4 1];        % 顶点2到顶点3、顶点9的边的权重
-a(3,4)=1;       % 同上。因为写过1到3，和2到3的边的权重，无需重复设
+a(3,4)=1;       % 同上。因为写过1到3，和2到3的边的权重，无需重复设
-a(4,5)=1;       
+a(4,5)=1;       
-a(5,6)=5;
+a(5,6)=5;
-a(6,7)=2; 
+a(6,7)=2; 
-a(7,8)=3; 
+a(7,8)=3; 
-a(8,9)=5;
+a(8,9)=5;
-s=cellstr(strcat('节点',int2str([1:9]')));
+s=cellstr(strcat('节点',int2str([1:9]')));
-G=graph(a,s,'upper');
+G=graph(a,s,'upper');
-P=plot(G,'layout','force');
+P=plot(G,'layout','force');
-T=minspantree(G,'Method','sparse');
+T=minspantree(G,'Method','sparse');
-L=sum(T.Edges.Weight);
+L=sum(T.Edges.Weight);
-highlight(P,T,"EdgeColor","red",'LineWidth',2.5);
+highlight(P,T,"EdgeColor","red",'LineWidth',2.5);

盘荣博/README.md

@@ -0,0 +1,7 @@
+**三个文件夹**
+
+数据可视化：暂时没有文件。
+
+图：其中包括最短路径和最小生成树的展示
+
+线性规划: 使用scipy库，便且使用scipy的optimize模块下的linprog函数，numpy函数统一操作数据以及使用matplotlib库下的pyplot模块

盘荣博/图/例10.py

@@ -0,0 +1,28 @@
+import networkx as nx
+import pylab as plt
+import numpy as np
+p=[25,26,28,31] ;a=[10,14,18,26] ; r=[20,16,13,11]
+b= np.zeros((5,5))
+
+for i in range(5):
+    for j in range(i+1,5):
+        b[i,j]=p[i]+np.sum(a[0:j-i])-r[j-i-1]
+G=nx.DiGraph(b)
+print(G)
+p=nx.dijkstra_path(G,source=0,target=4,weight='weight')
+print("最短路径为：",np.array(p)+1)#下标从零开始
+d=nx.dijkstra_path_length(G,source=0,target=4,weight="weight")
+print("所求的费用最小值为：",d)
+s=dict(zip(range(5),range(1,6)))
+plt.rc("font",size=16)
+pos=nx.shell_layout((G))#设置布局
+print(type(pos),'\npos=',pos)
+w=nx.get_edge_attributes(G,"weight")
+print(type(w),'\nw=',w)
+nx.draw(G,pos,font_weight='bold',labels=s,node_color='r')#绘制点和边
+nx.draw_networkx_edge_labels(G,pos,edge_labels=w)#绘制标签
+#绘制最短路径
+path_edges=list(zip(p,p[1:]))
+print(type(path_edges),"\npath_edges=",path_edges)
+nx.draw_networkx_edges(G,pos,edgelist=path_edges,edge_color="r",width=1)
+plt.savefig("figure10_9.png",dpi=1000);plt.show()

盘荣博/数据可视化/week1.py

@@ -0,0 +1,70 @@
+import pandas as pd
+<<<<<<< HEAD
+import numpy as np
+from  scipy.stats import zscore
+from sklearn.decomposition import PCA
+=======
+from  scipy.stats import zscore
+>>>>>>> remotes/origin/盘荣博
+import matplotlib.pyplot as plt
+from matplotlib.pyplot import ylabel
+df = pd.read_excel("棉花产量论文作业的数据.xlsx")
+# plt.plot(df["年份"],df["单产"])
+plt.rcParams['font.sans-serif']="SimHei"
+# plt.rcParams['size'] =10
+# plt.ylabel('单产')
+# plt.xlabel('年份')
+
+# print(df)
+d = df.to_numpy()[:,1:]
+print(d)
+plt.subplot(4,1,1)
+plt.scatter(d[:,:1],d[:,1:2],c='r')
+ylabel('原始数据'),plt.title("单产和种子费用的关系")
+#公式调用标准化，遵守标准正态分布
+data = zscore(d)
+print(data)
+plt.subplot(4,1,2)
+plt.scatter(data[:,:1],data[:,1:2],c='b',)
+ylabel('zscore')
+
+print(d.max(axis=0))
+print(d.std(axis=0))
+print(d.mean(axis=0))
+#手写标准正态分布
+data1=(d-d.mean(axis=0))/d.std(axis=0)
+print(data1)
+plt.subplot(4,1,3)
+plt.scatter(data1[:,:1],data1[:,1:2],c='y')
+ylabel('手写标准正态分布')
+
+data2=(d-d.min(axis=0))/(d.max(axis=0)-d.min(axis=0))
+plt.subplot(4,1,4)
+plt.scatter(data2[:,:1],data2[:,1:2],c='g')
+plt.xlabel('压缩到0~1')
+print(data==data1)
+
+<<<<<<< HEAD
+# plt.savefig("shuju.jpg",dpi=2000)
+# plt.show()
+md= PCA().fit(data)
+cf = np.cov(data.T)#求协方差矩阵
+print(cf)
+c, d= np.linalg.eig(cf)
+print("特征值：\n",c)
+print(md.explained_variance_)
+e=c/c.sum()
+# for _ in range(len(e)):
+#     if(_!=0):
+#         e[_]+=e[_-1]
+print('贡献率：')
+print(e)
+print(md.explained_variance_ratio_)
+print('特征向量：')
+print(d.T)
+print(md.components_)
+print(md.components_-d.T<=0.1)
+=======
+plt.savefig("shuju.jpg",dpi=2000)
+plt.show()
+>>>>>>> remotes/origin/盘荣博

盘荣博/线性规划/optimiaze模块/5.4投资的收益和风险模型1.py

@@ -0,0 +1,19 @@
+from matplotlib import pyplot  as plt
+from numpy import ones , diag , c_ , zeros
+from scipy.optimize import  linprog
+import time
+start = time.time()
+c=list( [-0.05,-0.27,-0.19,-0.185,-0.185])
+A=c_[zeros(4),diag([0.025,0.015,0.055,0.026])]
+Aeq = [[1,1.01,1.02,1.045,1.065]];beq=[[1]]
+a=0;aa=[];ss = []
+while a<=0.05:
+    b=list(ones(4)*a)
+    res= linprog(c,A,b,Aeq,beq)
+    aa.append(a);ss.append(-res.fun)
+    a=a+0.001
+end = time.time()
+print("花费时间：",end-start)
+plt.plot(aa,ss,"r*")
+plt.xlabel("$a$");plt.ylabel("$Q$",rotation=90)
+plt.savefig("figure5_1_1.png",dpi=500) ;plt.show()

盘荣博/线性规划/optimiaze模块/5.4投资的收益和风险模型2.py

@@ -0,0 +1,26 @@
+import numpy as np
+from numpy import ones , zeros , c_,diag
+from scipy.optimize import linprog
+import matplotlib.pyplot as plt
+c = np.append(zeros(5).tolist(),[1]).tolist()
+print(c)
+A=np.append(zeros(4).reshape(4,1),diag([0.025,0.015,0.055,0.026]),axis=1)
+A=np.append(A,ones(4).reshape(4,1)*-1,axis=1).tolist()
+Aeq =[[1,1.01,1.02,1.045,1.065,0]] ;beq=[1]
+A.append([-0.05,-0.27,-0.19,0.185,-0.185,0])
+print(A)
+k=0.05;step = 0.005
+b=([0]*4);b.append(-k)
+print(b)
+kk=[];ss=[]
+while k<0.28:
+    res= linprog(c,A,b,Aeq,beq)
+    kk.append(k)
+    ss.append(res.fun)
+    print(res.fun)
+    k+=step
+    b[4]=-k
+plt.plot(kk,ss,'r*')
+plt.xlabel("$k$");plt.ylabel('$R$')
+plt.savefig("figures5_1_2.png",dpi=500);plt.show()
+

盘荣博/线性规划/optimiaze模块/scipy函数.py

@@ -0,0 +1,56 @@
+from scipy.optimize import linprog
+import time
+import numpy as np
+'''
+目标函数：
+min z= -x1 + 4 x2 
+
+约束条件：
+-3x1 + x2 < = 6
+x1 + 2x2 <= 4 
+x2 >= -3 
+
+'''
+begin = time.time()
+c=[-1,4]
+A=[[-3,1],[1,2]]
+b=[[6],[4]]
+bounds=((None , None),(-3,None))
+res= linprog(c,A,b,None,None,bounds)
+end= time.time()
+print("结果为：",res.fun)
+print("x的值位",res.x," ",type(res.x))
+print(end-begin,"\n",'----------------------------------------')
+'''
+线性规划问题
+ max z = x1  - 2x2 - 3x3 
+ 约束条件
+ -2x1 + x2 + x3 <= 9
+ -3x1 + x2 + 2x3 >= 4 
+ 4x1 -2x2 -x3 = -6 
+ x1 >= -10 ,x2 >=0
+'''
+'''
+scipy标准形式：
+min w= -x1 + 2x2  + 3x3 
+
+约束： 
+-2x1 + x2 + x3 <= 9
+3x1 - x2 - 2x3 <= -4
+'''
+c= (-np.array([1,-2,-3])).tolist()
+A=[[-2,1,1],[-3,1,2]]
+A[1]= (-np.array(A)[1]).tolist()
+b= [[9],[4]]
+b[1]=(-np.array(b[1])).tolist()
+Aeq=[[4,-2,-1]]
+beq =[[-6]]
+LB =[-10,0,None]
+UB = [None]*len(c)
+bounds= tuple(zip(LB,UB))
+res= linprog(c,A,b,Aeq,beq,bounds)
+end= time.time()
+print("结果为：",res.fun)
+print("x的值位",res.x," ",type(res.x))
+print(end-begin,"\n",'----------------------------------------')
+print(res)

盘荣博/线性规划/optimiaze模块/例 5.3.py

@@ -0,0 +1,23 @@
+import time
+
+import numpy as np
+from scipy.optimize import linprog
+start = time.time()
+c = [-110,-120,-130,-110,-115,150]
+c = (-np.array(c)).tolist()
+A=[[1,1,0,0,0,0],
+   [0,0,1,1,1,0],
+   [8.8,6.1,2.0,4.2,5.0,-6],
+   [8.8,6.1,2.0,4.2,5.0,-3]
+   ]
+A[3]=(-np.array(A[3])).tolist()
+b=[[200],[250],[0],[0]]
+Aeq =[[1,1,1,1,1,-1]]
+beq =  [[0]]
+LB= [0]*len(c)
+UB= [None]*len(c)
+bounds= tuple(zip(LB,UB))
+res = linprog(c,A,b,Aeq,beq,bounds)
+end = time.time()
+print("最优解：\n",res.x)
+print("目标函数最小值：\n",-res.fun)

盘荣博/线性规划/插值和拟合/inter.py

@@ -0,0 +1,16 @@
+import numpy as np
+import pylab as pl
+from scipy import interpolate
+import matplotlib.pyplot as plt
+x = np.linspace(0,2*np.pi+np.pi/4,10)
+x1 = np.linspace(0,2*np.pi+np.pi/4,100)#num个0~2*Pi+Pi/4的范围的点
+y = np.sin(x)
+y1= np.sin(x1)
+plt.xlabel(f"安培/A")
+plt.ylabel(f'伏特/V')
+linear_ = interpolate.interp1d(x,y)
+print(linear_(x1))
+plt.rcParams['font.sans-serif'] = ['SimSun']#设置字体
+plt.plot(x,y,'o',label=f"原始数据")
+plt.plot(x1,linear_(x1),"*",label="线性插入")
+pl.show()

盘荣博/线性规划/插值和拟合/test1.py

@@ -0,0 +1,27 @@
+from scipy import  optimize
+import numpy as np
+
+c= np.array([-1,4])
+A=np.array([[-3,1],[1,2]])
+b=np.array([6,4])#小于关系
+Aeq=np.array([[1,1,1]])#相等
+beq = np.array([7])#
+bounds = ((None,None),(-3,None))
+res = optimize.linprog(c,A,b,None,None,bounds=bounds)
+print("目标函数最小值：",res.fun)#目标函数最优解
+print("最优解",res.x)#求得的最优解
+
+
+c = [-1,2,3]
+A=[[-2,1,1],[3,-1,-2]]
+b=[[9],[-4]]
+Aeq =[[4,-2,-1]]
+beq=[-6]
+LB=[-10,0,None]
+UB = [None]*len(c)
+print(UB)
+bound = tuple(zip(LB,UB))
+print(zip(LB,UB),'\n',bound)
+res = optimize.linprog(c,A,b,Aeq,beq,bounds=bound)
+print("函数的最小值为",res.fun)
+print("最优解为：",res.x)