: m/ y& d) |% K( j' n5 q 图3 MATLAB绘图面板中的图例 9 N6 H& k$ ^' a) H& S: J) s9 Z- K; o
要注意,需要在工作区选中变量后绘图面板中的这些图标才会激活。接下来就可以选中一个中意的图标进行绘图,一般都直接先选第一个(plot)看一下效果,然后再浏览整个面板,看看有没有更合适的。下面我们进行绘图操作。 . Q" e- u: W. C' y8 e2 c/ H / C G8 n5 Q; O& SStep2.2:选中变量 DataNum 和 Pclose,在绘图面板中单机 plot 图标,马上可以得到这两个变量的可视化结果,如图 4 所示,同时还可以在命令窗口区看到绘制此图的命令: 3 O: O) k+ e" T4 [/ I6 Z5 V ) A7 D6 a0 w) o/ u+ x9 B>> plot(DateNum,Pclose)8 X5 a! X \' f" t; G' @; O
7 K5 ~/ L9 D8 n2 i1 Q
4 J3 _4 m" C4 W8 g1 P : G+ I$ R" e& Q( k1 o- [ 图4 通过 plot 图标绘制的原图 ( p Q2 A/ z0 m# o& a8 x4 A8 [& A, V5 N0 N4 H" j0 _+ q
这样我们就知道了,下次再绘制这样的图直接用 plot 命令就可以了。一般情况下,用这种方式绘图的图往往不能满足我们的要求,比如我们希望更改: D$ N; S7 [: f. h l) s
; ~$ T! n7 b9 [# ]! }
(1)曲线的颜色、线宽、形状; & ?1 C6 ~9 X f( S- g) h + B2 v6 r H- e(2)坐标轴的线宽、坐标,增加坐标轴描述;6 k7 Y* L7 J. V+ K
* ?/ O- G% G% k( g
(3)在同个坐标轴中绘制多条曲线。 & x7 V Y2 F8 f6 F% d5 I2 E Z0 t z0 a$ E
此时我们就需要了解更多关于命令 plot 的用法,这时就可以通过 MATLAB 强大的帮助系统来帮助我们实现期望的结果。最直接获取帮助的两个命令是 doc 和 help,对于新手来说,推荐使用 doc,因为 doc 直接打开的是帮助系统中的某个命令的用法说明,不仅全,而且有应用实例,这样就可以“照猫画虎”,直接参考实例,从而将实例快速转化成自己需要的代码。 " ~* m- D. ]0 t& c$ M - S, ]: f! w5 q, y. o接下来我们就要考虑如何评估股票的价值和风险呢?9 @0 e( F# X/ V( H. H- Q
+ v" s! j- u2 y0 L/ v& E z
对于一只好的股票,我们希望股票的增幅越大越好,体现在数学上,就是曲线的斜率越大越好。; ?& X5 @) r) k) p; n
5 E6 v6 i2 K9 N, I# C 对于风险,则可用最大回撤率来描述更合适,什么是最大回撤率?! O5 U+ s! P1 L8 L# {, ]
; x" W" \* N- s0 A( Y
最大回撤率的公式可以这样表达:. P, B1 r. t( K3 r) E, a' q* J
: w3 B, c6 D: H/ Q% ]5 U- X1 `4 Z
D为某一天的净值,i为某一天,j为i后的某一天,Di为第i天的产品净值,Dj则是Di后面某一天的净值 ) [1 _6 G4 u/ q4 T8 {6 M& Y ) w ~! V w% ~2 `$ t, ndrawdown=max(Di-Dj)/Di,drawdown就是最大回撤率。其实就是对每一个净值进行回撤率求值,然后找出最大的。可以使用程序实现。最大回撤率越大,说明该股票的风险越高。所以最大回撤率越小,股票越好。% x7 H5 d& n' t, H( z
7 {. C4 x+ U' w! ?
斜率和最大回撤率不妨一个一个来解决。我们先来看如何计算曲线的斜率。对于这个问题,比较简单,由于从数据的可视化结果来看,数据近似成线性,所以不妨用多项式拟合的方法来拟合该改组数据的方程,这样我们就可以得到斜率。 ; f0 q( I9 a2 ?5 ~: O( J q$ B 9 b* q! ]- ?" Y* N, l9 J @7 @' IStep2.3:通过polyfit()多项式拟合的命令,并计算股票的价值,具体代码为: 0 b& p7 w8 O; N' [" q- P' c/ g" k! R
>> p = polyfit(DateNum,Pclose,1); % 多项式拟合 4 n7 d& V' v& U! W; a5 K1 O& a 5 |# w9 m$ k5 ]* g& [. f& E" G>> value = p(1) % 将斜率赋值给value,作为股票的价值 " n, E) G7 Q4 s( x/ G5 x! G4 P1 e5 n+ l* W3 A& b/ w' D
value =* z" w8 G- o0 z) R7 j
* A& q8 X% }9 b- r. e 0.12128 i8 J" {- L) m' ~4 U
0 h& J8 y. _0 j: a
代码分析:%后面的内容是注释。polyfit()有三个参数,前两个大家都能明白是什么意思,那第三个参数是什么意思呢?它表示多项式的阶数,也就是最高次数。比如:在本例中,第三个参数为1,说明其为一次项,即一次函数。第三个参数为你要拟合的阶数,一阶直线拟合,二阶抛物线拟合,并非阶次越高越好,看拟合情况而定。polyfit()返回阶数为 n 的多项式 p(x) 的系数,p 中的系数按降幂排列。在本例中的P(1)指的是最高项的系数,即斜率。) Q# v% @ g2 f. r8 P, A2 h, K* o
M" g, {, \% {7 c5 b8 r5 M3 c) [# X%% 输入数据0 e, b2 [; V- P7 x$ U
clc, clear all, close all 3 {7 t% n( b1 @x = [1.5, 4.5, 7.5,10.5,13.5,16.5,19.5,22.5,25.5];7 p2 `8 x/ m! o7 j1 Z; e8 \
y = [7.0,4.8,3.6,3.1,2.7,2.5,2.4,2.3,2.2];- W$ G' S$ Z. ?) W" ?
plot(x, y, '*', 'linewidth', 1) % 这里的linewidth指的是散点大小% {0 }- \0 R: V/ g1 ]
set(gca,'linewidth',2) % 设置坐标轴的线宽为2( T6 d9 d6 i* E! L v7 x4 t
xlabel('销售额x/万元','fontsize',12)) M& A ^2 \* m7 h( U
ylabel('流通率y/%','fontsize',12) 7 I$ i' ^+ Z1 ~ p(2)对数形式非线性回归0 r. T( K0 a A
# L4 m4 U6 ^- A3 G
%% 对数形式非线性回归. B w' Z" _# j7 T( @* q9 B* N
m1 = @(b,x) b(1) + b(2)*log(x);. ]6 U- c" W( x1 k1 x) o
nonlinfit1 = fitnlm(x,y,m1,[0.01;0.01]) 1 ^2 P$ P" E- `- C- kb = nonlinfit1.Coefficients.Estimate; : @ r& n: _6 y8 T5 AY1 = b(1,1) + b(2,1)*log(x);) X/ n: J* ?8 ?
hold on " z; q# l1 X/ h# B$ V: ?4 Yplot(x, Y1, '--k', 'linewidth',2) 1 i2 H, }" H" k7 q运行结果如下:, \$ v7 S$ a" `% {3 E7 n" }
: K8 }6 z$ K9 w& a
nonlinfit1 =/ ~9 J1 d% u: _# l7 h
( @4 I a2 T2 H/ K2 y' N2 ANonlinear regression model:8 g9 [7 p# P6 ?) C( b. l% B5 h
6 i, ~+ w. z5 b1 X. a4 m; w
y ~ b1 + b2*log(x)% x+ n0 y7 [. U6 }- A
3 R# k5 g; J: [7 h0 qEstimated Coefficients:) m7 D& u, S* t# D- j+ c
- ^! f5 G8 W' Y, Y$ \2 ?) Q, C Estimate SE tStat pValue % z( |8 K4 o1 G3 w7 q) p
$ a C1 S9 r- u& o
b1 7.3979 0.26667 27.742 2.0303e-089 }; i$ `$ F4 n5 D2 E
, ?1 _ P9 f% G
b2 -1.713 0.10724 -15.974 9.1465e-07 6 z% B' R. G, Y* [* U 7 H: x8 d( t3 A- C# s+ a" QR-Squared: 0.973, Adjusted R-Squared 0.969 ; l3 F1 m( [- K ) X4 ~6 ^) ~( V$ ]+ w5 JF-statistic vs. constant model: 255, p-value = 9.15e-07 7 M: Z- d r8 |* L9 e7 |- }* o8 }- ~2 w3 i& z0 y1 p3 F
(3)指数形式非线性回归 , C! l/ k; T. Y. X , P0 ?3 ^+ f3 R& j: O%% 指数形式非线性回归4 B( g- I& l6 y& T! Q
m2 = 'y ~ b1*x^b2'; * {& V' y* f$ c4 c2 ?0 n3 @' bnonlinfit2 = fitnlm(x,y,m2, [1;1]) 0 b, ?, T x7 B Wb1 = nonlinfit2.Coefficients.Estimate(1,1); $ \8 v6 `8 K! N( zb2 = nonlinfit2.Coefficients.Estimate(2,1)/ l4 C8 r. ~" K9 ~0 X7 @3 g
Y2 = b1*x.^b2;/ z& R; i2 M, {; f) D1 O
hold on; / @1 ^8 B L+ s' x8 Y. V. Splot(x,Y2,'r','linewidth',2)8 C* M, r4 F& u. i* r7 o7 @
legend('原始数据','a+b*lnx','a*x^b') % 图例 - m2 E% M2 P+ k运行结果如下: % n( n5 ^4 l% a9 P$ Q3 W! o , [: a4 A. F2 j, j- j5 i( P. jnonlinfit2 =& w) P, G0 p# u4 S7 Q. o- S6 D
+ I+ L; c3 i2 H7 u! V, FNonlinear regression model: ! X0 ~. B# z* p% H2 H ) Y7 H7 ?3 M! g1 u4 k y ~ b1*x^b2 ( j. P+ Y: |' q . i L, ?# g( Y( T9 ^1 i3 I; `9 y$ HEstimated Coefficients: 0 j. I1 T$ v+ J- A: T7 s% c9 @ $ e1 T4 B: A7 `7 p5 j Estimate SE tStat pValue * f* G2 d3 R! o- u! ] . }. B+ P4 Y! E0 q6 `; e b1 8.4112 0.19176 43.862 8.3606e-103 \4 Z# t P4 Z2 g: J7 T: e
' Q5 t8 a4 D5 S3 q. d$ d4 _, Z( F
b2 -0.41893 0.012382 -33.834 5.1061e-093 e1 n0 H/ D$ h" _# y: i9 E
% E9 v$ P/ P- [( i5 S
R-Squared: 0.993, Adjusted R-Squared 0.992% y) f! n. r9 }# S
: V' a! ~& O% [8 w1 o
F-statistic vs. zero model: 3.05e+03, p-value = 5.1e-11 - v! M6 c# ?% V+ {" Z 3 a# a* i% Y% e4 r8 r8 E* c在该案例中,选择两种函数形式进行非线性回归,从回归结果来看,对数形式的决定系数为 0.973 ,而指数形式的为 0.993 ,优于前者,所以可以认为指数形式的函数形式更符合 y 与 x 之间的关系,这样就可以确定他们之间的函数关系形式了。 ! q3 K8 l. y# i. w7 D% l6 p- K6 P$ x7 X2 I5 j t) k
2.多元回归 9 H% _( C1 H; l2 j6 N# t ) J v! ?: [8 |- g9 R# X1.多元线性回归/ o$ v/ n# `* x/ G
/ m# m5 d, A: d, u! S2)F 检验法:当 F > F1-α(m,n-m-1) ,即认为因变量 y 与自变量 x1,x2,...,xm 之间有显著的线性相关关系;否则认为因变量 y 与自变量 x1,x2,...,xm 之间线性相关关系不显著。本例 F=67.919 > F1-0.05( 3,20 ) = 3.10。& A5 r- {9 u! o) q% w5 ?$ c) v
0 E! S# M: h |2 S( o- H
3)p 值检验:若 p < α(α 为预定显著水平),则说明因变量 y 与自变量 x1,x2,...,xm之间显著地有线性相关关系。本例输出结果,p<0.0001,显然满足 p<α=0.05。1 b1 f C" y+ u2 K/ ]
+ j5 O$ P: R" Z
以上三种统计推断方法推断的结果是一致的,说明因变量 y 与自变量之间显著地有线性相关关系,所得线性回归模型可用。s^2 当然越小越好,这主要在模型改进时作为参考。 # ?, ~$ R3 |* _0 M X* \0 A- p 9 f* |" ~) m! r. {1 e" V5 \, Y$ B3. 逐步回归2 W; X4 P5 I ~6 d1 l
3 e; l0 p7 w9 p( v' n, F
[ 例4 ] (Hald,1960)Hald 数据是关于水泥生产的数据。某种水泥在凝固时放出的热量 Y(单位:卡/克)与水泥中 4 种化学成品所占的百分比有关:, N& Q' M$ h# Q1 x# w3 N/ d, k7 U" e
" y5 a0 R% z c# r0 e . v1 m: A, V. J% u8 r8 z9 ?) W : |# G U, ~: l! F/ w在生产中测得 12 组数据,见表5,试建立 Y 关于这些因子的“最优”回归方程。 - O) i, S* s2 r' q+ A: F: Q) r/ h; K5 ~
\8 a* i5 W4 E! Q# v |' `. @( n6 ^) ?; t0 I- S6 ?9 o5 ^% B, n- U
对于例 4 中的问题,可以使用多元线性回归、多元多项式回归,但也可以考虑使用逐步回归。从逐步回归的原理来看,逐步回归是以上两种回归方法的结合,可以自动使得方程的因子设置最合理。对于该问题,逐步回归的代码如下: 5 `" N) B1 L0 h' u, {; F1 b7 t# H) G% }8 ~" r$ J1 F
%% 逐步回归 , g! Q! w# L' m) N6 ]) dX=[7,26,6,60;1,29,15,52;11,56,8,20;11,31,8,47;7,52,6,33;11,55,9,22;3,71,17,6;1,31,22,44;2,54,18,22;21,47,4,26;1,40,23,34;11,66,9,12]; %自变量数据" ?: X) f1 [. Q% U
Y=[78.5,74.3,104.3,87.6,95.9,109.2,102.7,72.5,93.1,115.9,83.8,113.3]; %因变量数据$ ^& M( @* ?7 V# {; [% v
stepwise(X,Y,[1,2,3,4],0.05,0.10)% in=[1,2,3,4]表示X1、X2、X3、X4均保留在模型中 8 e* T) F/ P3 l( L* o! E) e" i程序执行后得到下列逐步回归的窗口,如图 4 所示。 V* \ |2 y' w( n9 Z9 E
: [. ]" C. I* s. }5 h4 N8 Q5 k6 e$ Q2 N9 Q' ]3 X1 Q7 C/ ?; ^
( [+ _, }" r, J
图4 4 v# g6 R. [: f; S, B9 K( w4 X* r6 Y% V
在图 4 中,用蓝色行显示变量 X1、X2、X3、X4 均保留在模型中,窗口的右侧按钮上方提示:将变量X4剔除回归方程(Move X4 out),单击 Next Step 按钮,即进行下一步运算,将第 4 列数据对应的变量 X4 剔除回归方程。单击 Next Step 按钮后,剔除的变量 X3 所对应的行用红色表示,同时又得到提示:将变量 X3 剔除回归方程(Move X3 out),单击 Next Step 按钮,这样一直重复操作,直到 “Next Step” 按钮变灰,表明逐步回归结束,此时得到的模型即为逐步回归最终的结果。最终结果如下: ! | P+ Y9 W, _4 s0 N# T, w. T2 Q& H$ E3 ~' M
. j& F% Y' R, m5 Y 0 j- B8 p3 _5 _1 A+ z+ U4. 逻辑回归 G. o* _/ q2 g R* G; a: `) C) E" ^4 z; J; Z5 Q, u9 Y1 g
[ 例5 ] 企业到金融商业机构贷款,金融商业机构需要对企业进行评估。评估结果为 0 , 1 两种形式,0 表示企业两年后破产,将拒绝贷款,而 1 表示企业 2 年后具备还款能力,可以贷款。在表 6 中,已知前 20 家企业的三项评价指标值和评估结果,试建立模型对其他 5 家企业(企业 21-25)进行评估。, n4 t& Z, P! C8 \8 Y
. F/ V+ o# T3 z( |5 B) Z : w. d6 a. c: H# `; p7 X. @' R8 u , r; h7 L+ v8 k$ n U对于该问题,很明显可以用 Logistic 模型来回归,具体求解程序如下:5 G- D( K7 P. ^" k! Z% I0 O. m
( W- i: Z/ x7 ]. Q$ ^& c# v
程序中需要用到的数据文件logistic_ex1.xlsx已上传github:https://github.com/xiexupang/mathematical-modeling/tree/master/%E5%9B%9E%E5%BD%92/%E9%80%BB%E8%BE%91%E5%9B%9E%E5%BD%92 ! R. {* u; @% j9 n# G& X % ]% Q! m4 i+ I: E7 G( \4 r% logistic回归4 F) o# R% y! H5 T2 \% _
6 g! Z9 L& }; N# d3 d
%% 导入数据 7 S' _4 |# B8 fclc,clear,close all . ?& S7 R% e) V8 G* g* Y" hX0 = xlsread('logistic_ex1.xlsx','A2:C21'); % 前20家企业的三项评价指标值,即回归模型的输入 3 R. r& K: l/ C( P2 yY0 = xlsread('logistic_ex1.xlsx','D221'); % 前20家企业的评估结果,即回归模型的输出: z; Q* g8 L a1 k/ X: [4 O' N& }
X1 = xlsread('logistic_ex1.xlsx','A2:C26'); % 预测数据输入 7 ^6 T$ ~1 W/ o2 N& Q+ o5 s, j* G9 n3 N
%% 逻辑函数1 t% P. H. I7 m$ ?$ I
GM = fitglm(X0,Y0,'Distribution','binomial');5 @! `2 n0 ?4 N( ^# D0 [* s
Y1 = predict(GM,X1); ! [3 H! G2 m# b- w* Y0 G& G ' ^+ P# K0 h* H7 e5 K- I8 e%% 模型的评估% q1 j- F$ V8 O( m8 d
N0 = 1:size(Y0,1); % N0 = [1,2,3,4,……,20] 5 v+ L( w% K5 m4 F2 KN1 = 1:size(Y1,1); % N1 = [1,2,3,4,……,25] 7 n6 | ^0 M) aplot(N0',Y0,'-kd'); % N0'指的是对N0'进行转置,N0'和Y0的形式相同,该行代码绘制的是前20家企业的评估结果6 u1 I2 I, T7 G/ {: q) w$ `
% plot()中的参数'-kd'的解析:-代表直线,k代表黑色,d代表菱形符号 $ a6 \3 z2 g7 ^) |- Hhold on;+ f8 w0 @$ h0 P8 C
scatter(N1',Y1,'b'); % N1'指的是对N1'进行转置,N1'和Y1的形式相同 1 M4 y6 z' M$ k6 Sxlabel('企业编号'); 5 C0 d) v( R$ ?9 R. r6 R- Gylabel('输出值');2 A0 H8 O1 E' K# M
得到的回归结果与原始数据的比较如图5所示。6 v3 G' ?. r- k- k7 }
; i, I g: y- L; ^6 O- `, c7 f
; r9 e+ H: c @( |8 ]