相关资料：
可以使用Matlab Online ，可以申请它的试用30天，不需要安装软件到本机，直接使用online 版本,但是online 版本的缺点是延迟高
推荐使用本地版（速度流畅，要求电脑配置高），可以自行寻找资源下载
Matlab 官方文档： https://www.mathworks.com/help/
课程视频：https://www.bilibili.com/video/BV1GJ41137UH?p=1
课件pdf : https://pan.baidu.com/s/1aT2B7B_mqyKMfJn-QolPiA 提取码: yswu
Github （个人代码）: https://github.com/zzzcb/matlab
注： 本博客是基于 Matlab 2015a 来的，所以不支持一些高版本的语法，例如后面高版本可以用双引号表示字符串，这里的2015a 只能用单引号。
注： 本博客中的参考代码是直接网上的答案，本人未经过验证，所以仅供参考！！！
Lecture 01 学习导引：
Lecture 02 基本操作与矩阵操作：
Introduction
基本界面：
有两种编程方式：
MATLAB as calculator
^:这里是次方号 ！
Note:
Command Window 如何清屏： clc 命令
遇到不熟悉的命令，可以使用谷歌，可以使用Matlab 自带文档，也可以在命令行中使用help ,
例如：在Command Window中：
e 如何表示: exp(1) ， e 的平方是 exp(2)
总结起来就是：遇到不会就online help ,
变量：
变量类型：
Matlab 中的关键字：
变量命名不规范容易导致下面的问题：
Note : clear + 变量名 会将工作区的该变量给清除， 如果只有clear 会清空工作区的所有变量 。
format 表示
>> format rat
>> 3/13 + 4/14 + 5/15
ans =
232/273
>> format long
>> 3/13+4/14 + 5/15
ans =
0.849816849816850
分号的作用：
一些常用的 指令：
Array operation
Array 可以分为 矩阵（Matrix）和向量（Vector）。
如何输入A矩阵呢？
数组索引：
大概是两种方法：
1，(row ，column) 用逗号分开
2，(idx) 按列索引
小练习：
答案： A(1,2) = 76 ; A(3,2) = 76; A([1 2],[2 3]) = 0;
冒号操作：
练习：如何去掉下面矩阵的第三行？
1, A= A([1 2],:)
2, A(3,:) = []
数组拼接：
三种：
1，空格 或 逗号 横着拼接
2，分号 竖着拼接
A =
1 2
3 4
B =
5 6
7 8
>> C = [A B] // 空格
C =
1 2 5 6
3 4 7 8
>> C = [A,B] // 逗号
C =
1 2 5 6
3 4 7 8
>> C = [A;B] // 分号 相当于换行
C =
1 2
3 4
5 6
7 8
答案：
A =
1 2
3 4
B =
9 9
9 9
C =
5 6 7 8
D =
-2 -1 0 1
>> F = [A B;C;D]
F =
1 2 9 9
3 4 9 9
5 6 7 8
-2 -1 0 1
数组操作：
主要有
+ 加
- 减
* 乘
/ 除 ，要用逆矩阵
^ 次方
. 点乘 或 点左除 或点右除 或点次方 dot
' 转置
一些特殊的矩阵：
注：
diag(A) 得到 矩阵A的对角向量。
rand() 返回一个 （0,1）的小数
与矩阵相关的常用函数：
注：
find(A>2) ,find 的参数一般跟一个 条件表达式，
Lecture 03 结构化程序与自定义函数：
Script writing
Structured programming
注： matlab 中的 ！= 是 ~=
if 语句：
switch 语句：
while 循环语句：
for 循环
预先申请空间：
预申请空间比动态内存更有效率。
break 和 continue
有用的小技巧
Functions
Note: which Locate functions and files. edit Edit or create a file
User-defined function
函数名一般和文件名一致，一般一个函数自己一个脚本，然后 其他脚本里面可以调用它，
有多个返回值：
函数内置变量：
等等，其他自行探索...
匿名函数：
总结：随用随定义，只能在 定义匿名函数的 文件中使用，
函数指针：
%% 函数 做为一个参数 f1 = @(x) sin(x); f2 = @(f,arg) f(arg); % 参数 f 为一个函数指针 ret = f2(f1,pi/6); disp(ret); ret2 = f2(@cos,pi/6); % @cos 引用内置的函数 disp(ret2);
input ，num2str ,isempty 的简单用法：
Lecture 04 数据结构和文件访问：
Variables:string ,structure,cell
字符串
字符串的拼接
1 space, 2 ； 如下：
反转字符串：
结构体：
结构体数组：
结构体相关的函数
嵌套的结构体：
Cell :
另一种形式的结构体，
访问cell 的内容
使用（ ） 访问的是 pointer ,
使用 {} 访问的才是 具体的内容，
cell 相关的函数
num2cell Convert numeric array into cell array.
mat2cell Break matrix up into a cell array of matrices.
关于mat2cell :
多维数组
使用cat() 函数可以拼接两个数组
其实，以上的操作也可以用cat 来完成， 如下：
所以可以用 cat 来构建 多维数组，
reshaper() 函数
判断数据的类型
Data access
load() 和 save()
将 workspace 工作区的数据保存到文件中：
如何只存储 指定的变量呢？
可以参数指定要存储的变量名， 这里是存储了a 变量
将 文件中的数据加载到workspace
xlsread() 和 xlswrite()
xlsread()
它会自动去掉Head的部分，只保留数字，
那么如何得到 文字部分呢？
得到的head 是个cell
xlswrite()
底层 文件操作：
向文件写入：
Note： type 命令，List program file 。 相当于Linux的cat 命令， 把文件的内容显示出来。
从文件中读出：
这里使用的是 fscanf ,
最终得到的 res 为：
Lecture 05 基础绘图：
Basic plotting
plot() 函数：
x = 1:10;
y = sin(x);
plot(x,y);
plot(x); % a line. plot(x) will replace plot(x,y)!
%%
x2 = 1:10;
y2 = sin(x2);
plot(x2,y2);
hold on % change mode, add plot ,not replace. default value is replace!
plot(x2); % a line. plot(x) will be added!
plot(x2+y2);
注：这里展示的仅仅只是一部分，更多请查看文档。
x2 = 1:10; y2 = sin(x2); plot(x2,y2,'.-r'); % marker is . line is - color is red
legend() 函数：
默认的legend 是加在 右上角的，当然也可以改变它的位置，
legend(...,'Location',LOC) adds a legend in the specified location, LOC, with respect to the axes. LOC may be either a 1x4 position vector or one of the following strings: 'North' inside plot box near top 'South' inside bottom 'East' inside right 'West' inside left 'NorthEast' inside top right (default for 2-D plots) 'NorthWest' inside top left 'SouthEast' inside bottom right 'SouthWest' inside bottom left 'NorthOutside' outside plot box near top 'SouthOutside' outside bottom 'EastOutside' outside right 'WestOutside' outside left 'NorthEastOutside' outside top right (default for 3-D plots) 'NorthWestOutside' outside top left 'SouthEastOutside' outside bottom right 'SouthWestOutside' outside bottom left 'Best' least conflict with data in plot 这个不错哦 'BestOutside' least unused space outside plot
title() 函数 和 ?label：
? 可以是 x,y ,z
x3 = 1:0.3:10;
y3 = sin(x3);
hold on % add ,not replace
plot(x3,y3,'--r');
plot(x3,x3,'--g');
legend('y=sin(x)','y=x') % add legend
title('title example.Latex:e^{x}. t = 0 to 2\pi')% title , latex
xlabel('x label test,t = 0 to 10.test:e^x 2\pi')
ylabel('y label test')
关于如何获取latex 文本 ： https://latex.91maths.com/
text() and annotation() 函数
x3 = 1:0.3:10;
y3 = sin(x3);
hold on % add ,not replace
plot(x3,y3,'--r');
plot(x3,x3,'--g');
line([1 10],[0 0]) % line from (1,0) to (10,0)
line([4 4], [0 sin(4)]) % line from (4,0) to (4 ,sin(4))
line([5 5], [0 sin(5)]) % line from (5,0) to (5 ,sin(5))
legend('y=sin(x)','y=x') % add legend
title('Just a test.')% title , latex
xlabel('x label test')
ylabel('y label test')
t = '$$ \int_{4} ^{5} \sin(x) dx $$'; % latex format
text(1.5,5,t,'Interpreter','latex') % pos is (1.5,5)
annotation('arrow',[0.25 0.45],[0.55 0.2]) % creates an annotation object
%head pos: (0.25,0.55) .* (width,height)% 这里的坐标是 归一化坐标！
%tail pos: (0.45,0.2) .* (width,height)
当然，annotation 也可以画其他的 形状，具体看帮助文档，
'rectangle'
'ellipse'
'textbox'
'line'
'arrow'
'doublearrow' = two headed arrow
'textarrow' = arrow with text at tail end
Graphical object properties
Figure Adjustment
在知道如何调成第二个图之前，我们首先要先知道一些知识：
1，Graphical Objects
而且它们还有一定的层次关系，
2 Figure Properties
如何进入这里： 在图片界面，View -> Property Editor
改变某个Object的属性：
第一步，在创建的时候保存，例如 h = plot(x,y);
gca 和 gcf :
x4 = 1:0.3:10; y4 = sin(x4); h = plot(x4,y4); ret1 = gca; % current axes object ret2 = gcf; % current figure object disp(ret1); disp(ret2);
output:
allchild() 函数：
x4 = 1:0.3:10; y4 = sin(x4); h = plot(x4,y4); ret1 = allchild(gca); ret2 = allchild(gcf); disp(ret1); disp(ret2);
OUTPUT：
Line with properties:
Color: [0 0.4470 0.7410]
LineStyle: '-'
LineWidth: 0.5000
Marker: 'none'
MarkerSize: 6
MarkerFaceColor: 'none'
XData: [1x31 double]
YData: [1x31 double]
ZData: [1x0 double]
Show all properties
=========================
10x1 graphics array:
Menu (figMenuHelp)
Menu (figMenuWindow)
Menu (figMenuDesktop)
Menu (figMenuTools)
Menu (figMenuInsert)
Menu (figMenuView)
Menu (figMenuEdit)
Menu (figMenuFile)
Toolbar (FigureToolBar)
Axes
第二步，使用 get() 和 set() 进行操作：
使用 get
x4 = 1:0.3:10;
y4 = sin(x4);
h = plot(x4,y4);
disp('----line object-----')
get(h) % line object's properties
disp('----axes object-----')
get(gca) % axes object's properties
disp('----figure object-----')
get(gcf) % figure object's properties
----line object-----
AlignVertexCenters: 'off'
Annotation: [1x1 matlab.graphics.eventdata.Annotation]
BeingDeleted: 'off'
BusyAction: 'queue'
ButtonDownFcn: ''
Children: []
Clipping: 'on'
Color: [0 0.4470 0.7410]
CreateFcn: ''
DeleteFcn: ''
DisplayName: ''
HandleVisibility: 'on'
HitTest: 'on'
Interruptible: 'on'
LineStyle: '-'
LineWidth: 0.5000
Marker: 'none'
MarkerEdgeColor: 'auto'
MarkerFaceColor: 'none'
MarkerSize: 6
Parent: [1x1 Axes]
PickableParts: 'visible'
Selected: 'off'
SelectionHighlight: 'on'
Tag: ''
Type: 'line'
UIContextMenu: []
UserData: []
Visible: 'on'
XData: [1x31 double]
XDataMode: 'manual'
XDataSource: ''
YData: [1x31 double]
YDataSource: ''
ZData: [1x0 double]
ZDataSource: ''
----axes object-----
ALim: [0 1]
ALimMode: 'auto'
ActivePositionProperty: 'outerposition'
AmbientLightColor: [1 1 1]
BeingDeleted: 'off'
Box: 'on'
BoxStyle: 'back'
BusyAction: 'queue'
ButtonDownFcn: ''
CLim: [0 1]
CLimMode: 'auto'
CameraPosition: [5.5000 0 17.3205]
CameraPositionMode: 'auto'
CameraTarget: [5.5000 0 0]
CameraTargetMode: 'auto'
CameraUpVector: [0 1 0]
CameraUpVectorMode: 'auto'
CameraViewAngle: 6.6086
CameraViewAngleMode: 'auto'
Children: [1x1 Line]
Clipping: 'on'
ClippingStyle: '3dbox'
Color: [1 1 1]
ColorOrder: [7x3 double]
ColorOrderIndex: 2
CreateFcn: ''
CurrentPoint: [2x3 double]
DataAspectRatio: [4.5000 1 1]
DataAspectRatioMode: 'auto'
DeleteFcn: ''
FontAngle: 'normal'
FontName: 'Helvetica'
FontSize: 10
FontSmoothing: 'on'
FontUnits: 'points'
FontWeight: 'normal'
GridAlpha: 0.1500
GridAlphaMode: 'auto'
GridColor: [0.1500 0.1500 0.1500]
GridColorMode: 'auto'
GridLineStyle: '-'
HandleVisibility: 'on'
HitTest: 'on'
Interruptible: 'on'
LabelFontSizeMultiplier: 1.1000
Layer: 'bottom'
LineStyleOrder: '-'
LineStyleOrderIndex: 1
LineWidth: 0.5000
MinorGridAlpha: 0.2500
MinorGridAlphaMode: 'auto'
MinorGridColor: [0.1000 0.1000 0.1000]
MinorGridColorMode: 'auto'
MinorGridLineStyle: ':'
NextPlot: 'replace'
OuterPosition: [0 0 1 1]
Parent: [1x1 Figure]
PickableParts: 'visible'
PlotBoxAspectRatio: [1 0.7903 0.7903]
PlotBoxAspectRatioMode: 'auto'
Position: [0.1300 0.1100 0.7750 0.8150]
Projection: 'orthographic'
Selected: 'off'
SelectionHighlight: 'on'
SortMethod: 'childorder'
Tag: ''
TickDir: 'in'
TickDirMode: 'auto'
TickLabelInterpreter: 'tex'
TickLength: [0.0100 0.0250]
TightInset: [0.0506 0.0532 0.0134 0.0202]
Title: [1x1 Text]
TitleFontSizeMultiplier: 1.1000
TitleFontWeight: 'bold'
Type: 'axes'
UIContextMenu: []
Units: 'normalized'
UserData: []
View: [0 90]
Visible: 'on'
XAxisLocation: 'bottom'
XColor: [0.1500 0.1500 0.1500]
XColorMode: 'auto'
XDir: 'normal'
XGrid: 'off'
XLabel: [1x1 Text]
XLim: [1 10]
XLimMode: 'auto'
XMinorGrid: 'off'
XMinorTick: 'off'
XScale: 'linear'
XTick: [1 2 3 4 5 6 7 8 9 10]
XTickLabel: {10x1 cell}
XTickLabelMode: 'auto'
XTickLabelRotation: 0
XTickMode: 'auto'
YAxisLocation: 'left'
YColor: [0.1500 0.1500 0.1500]
YColorMode: 'auto'
YDir: 'normal'
YGrid: 'off'
YLabel: [1x1 Text]
YLim: [-1 1]
YLimMode: 'auto'
YMinorGrid: 'off'
YMinorTick: 'off'
YScale: 'linear'
YTick: [1x11 double]
YTickLabel: {11x1 cell}
YTickLabelMode: 'auto'
YTickLabelRotation: 0
YTickMode: 'auto'
ZColor: [0.1500 0.1500 0.1500]
ZColorMode: 'auto'
ZDir: 'normal'
ZGrid: 'off'
ZLabel: [1x1 Text]
ZLim: [-1 1]
ZLimMode: 'auto'
ZMinorGrid: 'off'
ZMinorTick: 'off'
ZScale: 'linear'
ZTick: [-1 0 1]
ZTickLabel: ''
ZTickLabelMode: 'auto'
ZTickLabelRotation: 0
ZTickMode: 'auto'
----figure object-----
Alphamap: [1x64 double]
BeingDeleted: 'off'
BusyAction: 'queue'
ButtonDownFcn: ''
Children: [1x1 Axes]
Clipping: 'on'
CloseRequestFcn: 'closereq'
Color: [0.9400 0.9400 0.9400]
Colormap: [64x3 double]
CreateFcn: ''
CurrentAxes: [1x1 Axes]
CurrentCharacter: ''
CurrentObject: []
CurrentPoint: [-192 -204]
DeleteFcn: ''
DockControls: 'on'
FileName: ''
GraphicsSmoothing: 'on'
HandleVisibility: 'on'
IntegerHandle: 'on'
Interruptible: 'on'
InvertHardcopy: 'on'
KeyPressFcn: ''
KeyReleaseFcn: ''
MenuBar: 'figure'
Name: ''
NextPlot: 'add'
Number: 1
NumberTitle: 'on'
PaperOrientation: 'portrait'
PaperPosition: [0.2500 2.5000 8 6]
PaperPositionMode: 'manual'
PaperSize: [8.5000 11]
PaperType: 'usletter'
PaperUnits: 'inches'
Parent: [1x1 Root]
Pointer: 'arrow'
PointerShapeCData: [16x16 double]
PointerShapeHotSpot: [1 1]
Position: [793 383 560 420]
Renderer: 'opengl'
RendererMode: 'auto'
Resize: 'on'
SelectionType: 'normal'
SizeChangedFcn: ''
Tag: ''
ToolBar: 'auto'
Type: 'figure'
UIContextMenu: []
Units: 'pixels'
UserData: []
Visible: 'on'
WindowButtonDownFcn: ''
WindowButtonMotionFcn: ''
WindowButtonUpFcn: ''
WindowKeyPressFcn: ''
WindowKeyReleaseFcn: ''
WindowScrollWheelFcn: ''
WindowStyle: 'normal'
使用 set 修改，
x4 = 1:0.3:10; y4 = sin(x4); h = plot(x4,y4); 图像如下：
现在要修改它的属性，将 x 限制在 [pi,3 pi] y 限制在 [-1.5 ,1.5 ]。
x4 = 1:0.3:10;
y4 = sin(x4);
h = plot(x4,y4); % line object
set(gca,'xlim',[pi,3*pi]) % x in [pi,3pi] set(gca,'ylim',[-1.5,1.5]) % y in [-1.5,1.5] Just a test!
Setting Font and Tick of Axes
x4 = 1:0.3:10;
y4 = sin(x4);
h = plot(x4,y4); % line object
set(gca,'xlim',[pi,3*pi]) % x in [pi,3pi]
set(gca,'ylim',[-1.5,1.5]) % y in [-1.5,1.5] Just a test!
set(gca,'fontsize',15) % set tick' font size equal to 15px
x4 = 1:0.3:10;
y4 = sin(x4);
h = plot(x4,y4); % line object
set(gca,'xlim',[pi,3*pi]) % x in [pi,3pi]
set(gca,'ylim',[-1.5,1.5]) % y in [-1.5,1.5] Just a test!
set(gca,'fontsize',15) % set tick' font size equal to 15px
set(gca,'xtick',pi:pi/2:3*pi) % use pi on tick. set(gca,'xticklabel',180:90:3*180)
Setting Line Style
x5 = 1:0.3:10;
y5 = sin(x5);
h = plot(x5,y5); % line object
set(h,'linestyle','-.',... 'linewidth',7.0,... 'color','g')
delete(h); 将会删除图像， delete： Delete file or graphics object.
修改Maker 的 Style:
x5 = 1:10; y5 = sin(x5); h = plot(x5,y5,'-d');
x5 = 1:10;
y5 = sin(x5);
h = plot(x5,y5,'-d'); % 设置Marker 的前提是，figure 首先要有 Marker，这里指定Marker 为d,即diamond.
set(h,'markeredgecolor','r',... 'markerfacecolor','g',... 'markersize',10) % set marker style 效果如下：
小练习：
完成图1 到 图2 的转换：
%% figure01
%x = 1:0.02*pi:2;% can't get 2!
x = linspace(1,2,90); % can get 2! num is 90,
y1 = x.^2;
y2 = sin(2*pi*x);
hold on; % add not replace!
l1 = plot(x,y1,'k');% color is black
l2 = plot(x,y2,'ro');% marker shape is circle,color is red
legend([l1,l2],'t^2','sin(2\pi t)','location','northwest')% legend
title('Mini Assignment #1')
xlabel('Time (ms)')
ylabel('f(t)')
set(gca,'xtick',linspace(1,2,6));%[1,2] num is 6
%% figure02
%x = 1:0.02*pi:2;% can't get 2!
x = linspace(1,2,90); % can get 2! num is 90,
y1 = x.^2;
y2 = sin(2*pi*x);
hold on; % add not replace!
l1 = plot(x,y1,'k');% color is black
l2 = plot(x,y2,'ro');% marker shape is circle,color is red
legend([l1,l2],'t^2','sin(2\pi t)','location','northwest')% legend
title('Mini Assignment #1')
xlabel('Time (ms)')
ylabel('f(t)')
set(gca,'xtick',linspace(1,2,6));%[1,2] num is 6
% to figure2 set(gca,'fontsize',16) % setting axes's fontsize set(l1,'linewidth',4) % setting l1's linewidth to 4 set(l2,'color',[160,32,240]/255,... % change color to purple 'marker','.',... % change marker to point 'markersize',16)
Multiple Figures
Note : 需要注意的是，此时会有两个figure，如果使用gcf gca 仅仅指的是最后一个figure，前面的figure就不能访问了，
%% multipe figures
x = 1:0.1:10;
y1 = x.^2;
y2 = exp(x);
%====
figure('position',[100,200,400,400]);% figure 1
plot(x,y1);
%====
figure('position',[600,600,400,400]);% figure 2
plot(x,y2);
Small figures in One Figure
%% small figures in one figure
x = linspace(-1,1,50);
y = sqrt(1-x.^2);
x = [x,fliplr(x)];%fliplr Flip array in left/right direction.%reverse
y = [y,-fliplr(y)];
subplot(221)
plot(x,y);
axis normal
title('normal[default]')
subplot(222)
plot(x,y);
axis square %makes the current axis box square in size. x axis's length == y axis's length
title('square')
subplot(223)
plot(x,y);
axis equal % x axis's unit == y axis's unit
title('equal')
%% small figures in one figure
x = linspace(-1,1,50);
y = sqrt(1-x.^2);
x = [x,fliplr(x)];%fliplr Flip array in left/right direction.%reverse
y = [y,-fliplr(y)];
subplot(221)
plot(x,y);
axis normal
title('normal[default]')
grid on box off
subplot(222)
plot(x,y);
axis square %makes the current axis box square in size. x axis's length == y axis's length
title('square')
axis off
subplot(223)
plot(x,y);
axis equal % x axis's unit == y axis's unit
title('equal')
grid on
我的使用习惯：
grid on
box off
保存图片到文件
Vector 格式 是矢量图 放大不失真！
%% save to file
x = linspace(1,10,100);
y = sin(x);
plot(x,y);
saveas(gcf,'test01','png') % bitmap format
saveas(gcf,'test01','svg') % vector format
更强大的保存图片参考 print ！
lecture 06 高级绘图
先上一张图：
Advanced 2D plots
x = logspace(-1,1,100); % yi 10 wei di
y = x.^2;
subplot(221);
plot(x,y);
title('plot');
subplot(222);
semilogx(x,y);
title('semilogx');
subplot(223);
semilogy(x,y);
title('semilogy');
subplot(224);
loglog(x,y);
title('loglog');
plotyy 函数
Graph with y-tick labels on the left and right side
%% plotyy() left y label and right y label
x = 1:0.1:10;
y1 = sin(x);
y2 = 2*sin(2*x) + 1;
[ax,h1,h2] = plotyy(x,y1,x,y2);% the handles of the two axes
% AX(1) is the left axes and AX(2) is the right axes.
set(get(ax(1),'ylabel'),'string','Left Y-axis')
set(get(ax(2),'ylabel'),'string','Right Y-axis')
title('Labeling plotyy');
set(h1,'linestyle','--');
set(h2,'linestyle',':'); % 如下图
Histogram 直方图
看分布是如何的，与概率相关，它是看整体的情况
%% histogram
x = randn(1,1000); % Normally distributed ,x is 1x1000 martix
subplot(211);
hist(x,10); % bins is 10. 10 groups
title('bins is 10')
subplot(212);
hist(x,50); % bins is 50. 50 groups
title('bins is 50') % 效果如下图
Bar Charts 柱状图 条形图
%% bar charts
x = 1:3; % 3 months
y = randi([2,10],1,3); % random 2-10. shape is 1x3 matrix.
subplot(221)
bar(x,y);
set(gca,'xticklabel',{'jan','feb','mar'}); % change xtick to text
title('bar single')
disp(x);
disp(y);
disp('===');
subplot(222)
y2 = randi([2,10],1,3);
y3 = randi([2,10],1,3);
bar(x,[y;y2;y3]); % The length of X must match the number of rows of Y.
title('bar multiple')
disp(x);
disp([y;y2;y3]);
disp('===');
subplot(223)
y = 1:3; % y axis
z = randi([2,10],1,3); % z axis
bar3(y,z); % y and z
title('bar3 single')
disp(y);
disp(z);
disp('===');
subplot(224)
y = 1:3; % y axis
z = randi([2,10],1,3); % z axis
z2 = randi([2,10],1,3); % z axis
z3 = randi([2,10],1,3); % z axis
bar3(y,[z;z2;z3]); % The length of Y must match the number of rows of Z.
title('bar3 multiple')
disp(y);
disp([z;z2;z3]);
disp('===');
OUTPUT
===
1 2 3
2 2 4
===
1 2 3
2 2 4
2 5 5
6 10 4
===
1 2 3
5 9 5
===
1 2 3
3 6 9
5 8 7
6 6 10
===
Stacked and Horizontal Bar Charts
%% stacked and horizontal bar charts
x = 1:3; % 3 months
y = randi([2,10],1,3); % random 2-10. shape is 1x3 matrix.
y2 = randi([2,10],1,3);
y3 = randi([2,10],1,3);
subplot(131)
bar(x,[y;y2;y3],'stacked');
title('stacked bar')
subplot(132);
barh(x,[y;y2;y3]);
title('horizontal bar');
subplot(133);
barh(x,[y;y2;y3],'stacked');
title('stacked horizontal bar');
disp(x);
disp([y;y2;y3]);
OUTPUT
1 2 3
6 6 7
8 8 4
8 7 3
Pie Charts
%% Pie Charts x = [10 5 20 30]; subplot(221); pie(x); subplot(222); pie(x,[0,0,0,1]);% 将最后一份分割出来 subplot(223); pie(x,[1,1,1,1]);% 每一份都分割开来 subplot(224); pie3(x,[0,0,0,1]);% 3d 绘图
Polar Chart 极坐标图
%% Polar Chart theta = 0:0.1*pi:pi/2; r = 2*cos(theta); subplot(121) polar(theta,r); theta2 = 0:0.1*pi:pi; r2 = 2*cos(theta2); subplot(122) polar(theta2,r2);
更多图形：
Stairs and Stem Charts 阶梯图和杆图
%% Stairs and Stem Charts
x = linspace(0,4*pi,40);
y = sin(x);
subplot(121);
stairs(x,y); % 阶梯图
title('stairs');
subplot(122);
stem(x,y); % 杆图
title('stem');
Boxplot and Error Bar 箱图 和 误差图
箱型图 略.
%% Error Bar x = 0:pi/10:pi; y = sin(x); err = std(y)*ones(size(x)); errorbar(x,y,err);
fill() 函数
%% fill theta = [1:2:15]*pi/8; x = cos(theta); y = sin(theta); subplot(121) plot(x,y);axis square; subplot(122) fill(x,y,'r');axis square;
Color space
set(lineObj,'color',[160,32,240]/255);
imagesc 函数， Scale data and display as image
%% imagesc d = [0 2 4 3; 8 10 12 14; 16 18 20 22]; imagesc(d) colorbar % 显示 color bar
d = [0 2 4 3; 8 10 12 14; 16 18 20 22]; imagesc(d) colormap(hot); % hot cool gray colorbar
内置colormaps
%% built-in colormaps
x = randi([0,255],10,1);
y = randi([0,255],1,10);
d = mod(x*y,255);
% imagesc(d); % default colormap is PARULA
%colormap('jet');
%imagesc(d);
colormap('hsv');
imagesc(d);
这个要自定义 ColorMap ,
%% custom colormap
x = [1:10; 3:12; 5:14];
myColorMap = zeros(256,3);
myColorMap(:,2) = (0:255)/255;
colormap(myColorMap);
imagesc(x);
colorbar;
3D plots
2D vs. 3D
%% 2d vs 3d
x=0:0.1:2*pi;
subplot(121);
plot(x,sin(x));
grid on
xlabel('x-axis'); ylabel('y-axis');
subplot(122);
plot3(x,sin(x),zeros(size(x)));
grid on
xlabel('x-axis'); ylabel('y-axis'); zlabel('z-axis');
Note :3维空间坐标系之---左手坐标系和右手坐标系的区别：
右手坐标系的z轴指向内，左手指向外。
Matlab 使用的是右手坐标系。
3D surface 图：
%% surface plots
x = linspace(-3,3,20);
y = linspace(-4,4,20);
[X,Y] = meshgrid(x,y);
Z = X.*exp(-X.^2-Y.^2);
subplot(121);
mesh(X,Y,Z);%3-D mesh surface. 3d 网格表面
subplot(122);
surf(X,Y,Z);%3-D colored surface.3d 颜色表面
contour() 函数
contour 是轮廓的意思，
Projection of equal heights of 3D plot onto a 2D plane
%% contour x = linspace(-3,3,20); y = linspace(-4,4,20); [X,Y] = meshgrid(x,y); Z = X.*exp(-X.^2-Y.^2); subplot(121); surf(X,Y,Z);%3-D colored surface. axis square; subplot(122); contour(X,Y,Z); axis square;
More Contour Plots
%% more contour x = linspace(-3,3,20); y = linspace(-4,4,20); [X,Y] = meshgrid(x,y); Z = X.*exp(-X.^2-Y.^2); subplot(121); surf(X,Y,Z);%3-D colored surface. axis square; subplot(122); [C,h] = contour(X,Y,Z); clabel(C,h); % 标上数值 axis square
meshc() and surfc() 函数
Combination of surface/mesh and contours。
%% meshc() and surfc() x = linspace(-3,3,20); y = linspace(-4,4,20); [X,Y] = meshgrid(x,y); Z = X.*exp(-X.^2-Y.^2); subplot(121); meshc(X,Y,Z); axis square; subplot(122); surfc(X,Y,Z); axis square;
View Angle: view()
观察的角度，
Azimuth and Elevation:
Note: AZ = -37.5, EL = 30 is the default 3-D view.
AZ = 0, EL = 90 is directly overhead and the default 2-D view.
%% view angle :view()
% Note: AZ = -37.5, EL = 30 is the default 3-D view.
% AZ = 0, EL = 90 is directly overhead and the default 2-D view.
x = linspace(-3,3,20);
y = linspace(-4,4,20);
[X,Y] = meshgrid(x,y);
Z = X.*exp(-X.^2-Y.^2);
subplot(121);
surf(X,Y,Z);
axis square;
subplot(122);
surf(X,Y,Z);
axis square;
view(0,90);
light() 函数：
light Create light.
打光的角度
%% light
x = linspace(-3,3,20);
y = linspace(-4,4,20);
[X,Y] = meshgrid(x,y);
Z = X.*exp(-X.^2-Y.^2);
subplot(221);
surf(X,Y,Z);
axis square;
xlabel('x');
ylabel('y');
zlabel('z');
subplot(222);
surf(X,Y,Z);
axis square;
light('position',[0,0,0.5]);
hold on;
plot3(0,0,0.5,'ro');
subplot(223);
surf(X,Y,Z);
axis square;
light('position',[5,5,0.5]);
hold on;
plot3(5,5,0.5,'ro');
subplot(224);
surf(X,Y,Z);
axis square;
light('position',[-5,-5,0.5]);
hold on;
plot3(-5,-5,0.5,'ro');
图中 红点为 打光的位置！
patch() 函数
patch 是补丁的意思
A graphical object containing polygons
%% patch()
x = [0,1,1,0];
y = [0,0,1,1];
subplot(231);
patch(x,y,'red')
x2 = [2,5; 2,5; 8,8];
y2 = [4,0; 8,2; 4,0];
subplot(232);
patch(x2,y2,'green')
v = [0, 0; 1 ,0; 1, 1; 0, 1];
f = [1,2,3]; % 指定要上颜色的面(3 个点)
subplot(233);
patch('faces',f,'vertices',v,'facecolor','red')
v = [0 0 0; 1 0 0 ; 1 1 0; 0 1 0; 0.25 0.25 1; ...
0.75 0.25 1; 0.75 0.75 1; 0.25 0.75 1];
f = [1 2 3 4; 5 6 7 8; 1 2 6 5; 2 3 7 6; 3 4 8 7; 4 1 5 8]; %
subplot(234);
patch('faces',f,'vertices',v, 'facecolor', 'cyan');
view(3); %3-d view
v = [0 0 0; 1 0 0 ; 1 1 0; 0 1 0; 0.25 0.25 1; ...
0.75 0.25 1; 0.75 0.75 1; 0.25 0.75 1];
f = [1 2 3 4; 5 6 7 8; 2 3 7 6; 3 4 8 7; 4 1 5 8]; % 没有这个面 1 2 6 5;
subplot(235);
patch('faces',f,'vertices',v, 'facecolor', 'cyan');
view(3); %3-d view
lecture 07 GUI 界面
如何开启？
设置显示组件的名称：
有疑问自行谷歌，（会QT的人感觉到很舒服）
对齐工具：
修改文字内容：
.fig 对应的 .m 文件：
% --- Executes on button press in pushbutton1. function pushbutton1_Callback(hObject, eventdata, handles) % hObject handle to pushbutton1 (see GCBO) % eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA) x = 1:0.1:10; y = sin(x); plot(x,y); % --- Executes on button press in pushbutton2. function pushbutton2_Callback(hObject, eventdata, handles) % hObject handle to pushbutton2 (see GCBO) % eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA) x = 1:0.1:10; y = cos(x); plot(x,y); % --- Executes on button press in pushbutton3. function pushbutton3_Callback(hObject, eventdata, handles) % hObject handle to pushbutton3 (see GCBO) % eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA) x = 1:0.1:10; y = tan(x); plot(x,y);
% --- Executes on button press in pushbutton1. function pushbutton1_Callback(hObject, eventdata, handles) % hObject handle to pushbutton1 (see GCBO) % eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA) x = 1:0.1:10; y = sin(x); plot(handles.axes1,x,y); % 指定在 axes1 上作图 % --- Executes on button press in pushbutton2. function pushbutton2_Callback(hObject, eventdata, handles) % hObject handle to pushbutton2 (see GCBO) % eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA) x = 1:0.1:10; y = cos(x); axes(handles.axes4); % 指定在 axes4 上作图 plot(x,y);
% --- Executes on slider movement. function slider1_Callback(hObject, eventdata, handles) % hObject handle to slider1 (see GCBO) % eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA) % Hints: get(hObject,'Value') returns position of slider % get(hObject,'Min') and get(hObject,'Max') to determine range of slider val1 = get(hObject,'Value'); val2 = get(handles.slider2,'Value'); set(handles.res_text,'String',num2str(int16(val1+val2))); % key -value % --- Executes on slider movement. function slider2_Callback(hObject, eventdata, handles) % hObject handle to slider2 (see GCBO) % eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA) % Hints: get(hObject,'Value') returns position of slider % get(hObject,'Min') and get(hObject,'Max') to determine range of slider val1 = get(hObject,'Value'); val2 = get(handles.slider1,'Value'); set(handles.res_text,'String',num2str(int16(val1+val2))); % key -value
全局变量
编译GUI程序 为 可执行文件：
需要使用 deploytool 指令，
选择第一个，
Lecture 08 Image Processing 图像处理 I：
Introduction to digital image
Types of Digital Image
三种具体的Image
Read and show images
Read and Show An Image
clear;
close all; % 关闭所有的figure
I = imread('pout.tif');% 图片为matlab 内置的图片
imshow(I); % show
我们可以通过改变这些数值来修改图片，
例如如下程序：
clear;
close all; % 关闭所有的figure
I = imread('pout.tif');% 图片为matlab 内置的图片
subplot(121);
imshow(I); % show
% 对图像的数据进行操作
for i = 1:size(I,1)%size(I)[1]
for j = 1:size(I,2)%size(I)[2]
if mod(i,2) == 0 && mod(j,2) == 0
I(i,j) = 0; % 0 是黑色的
end
end
end
subplot(122);
imshow(I); % show
imageinfo（） 函数
它可以查看图像的信息
%% imageinfo()
%I = imread('pout.tif');% 图片为matlab 内置的图片
imageinfo('pout.tif'); % 获取文件的相关信息
从中我们也可以看出 pout.tif 是内置的图像，可以通过路径把它找出来 ，
imtool() 函数
%% imtool()
imtool('pout.tif'); % isplay image in the Image Tool.
Image arithmetic（Image 的四则运算）
Image Multiplication: immultiply()
%% immultiply() image 乘法
I = imread('rice.png');% 内置的图片
subplot(121);
imshow(I);
J = immultiply(I,1.5); % 其实是将像素值 都乘以1.5
subplot(122);
imshow(J);
Image Addition: imadd()
%% imadd() image 加法
I = imread('rice.png');% 内置的图片
J = imread('cameraman.tif');% 内置
K = imadd(I,J);
subplot(131);
imshow(I);
subplot(132);
imshow(J);
subplot(133);
imshow(K);
注意，相加的时候，这两个图片的 dim 要相同。
Image Histogram: imhist() （直方图）
isplay histogram of image data.
%% imhist() isplay histogram of image data.
I = imread('pout.tif');% gray image
subplot(121);
imshow(I);
subplot(122);
imhist(I);
注：一般都是展示 灰度图的 直方图，
Histogram Equalization: histeq() 直方图均衡化
Enhances the contrast of the image ，增强图像的对比度
那么为什么要做 图像均衡化，
把集中在小区间的像素值 变为 大区间的！
Geometric Transformation 几何变换
这和前面不一样的是，这里改变的是 像素的位置， 前面改变的像素值。
一遍是通过矩阵进行改变：
左乘
更方便的是使用已经提供的 内置函数，
例如做旋转：
Image Rotation: imrotate()
%% Image rotate
I = imread('rice.png');
subplot(121);
imshow(I);
J = imrotate(I, 35, 'bilinear'); % 双线性插值 参见： https://www.cnblogs.com/linkr/p/3630902.html
subplot(122);
imshow(J);
size(I)
size(J)
注：看起来两个图片的大小一样，其实是不一样的，第二个的要比第一个大一些，
练习：不使用 imrotate() 实现图片的旋转，要求使用旋转矩阵
jiaodu=45; %要旋转的角度，旋转方向为顺时针
img=imread('rice.png'); %这里v为原图像的高度，u为原图像的宽度
subplot(121)
imshow(img); %这里y为变换后图像的高度，x为变换后图像的宽度
[h w]=size(img);
theta=jiaodu/180*pi;%转换为弧度
rot=[cos(theta) -sin(theta) 0;
sin(theta) cos(theta) 0;
0 0 1];
pix1=[1 1 1]*rot; %变换后图像左上点的坐标
pix2=[1 w 1]*rot; %变换后图像右上点的坐标
pix3=[h 1 1]*rot; %变换后图像左下点的坐标
pix4=[h w 1]*rot; %变换后图像右下点的坐标
height=round(max([abs(pix1(1)-pix4(1))+0.5 abs(pix2(1)-pix3(1))+0.5])); %变换后图像的高度
width=round(max([abs(pix1(2)-pix4(2))+0.5 abs(pix2(2)-pix3(2))+0.5])); %变换后图像的宽度
imgn=zeros(height,width);
delta_y=abs(min([pix1(1) pix2(1) pix3(1) pix4(1)])); %取得y方向的负轴超出的偏移量
delta_x=abs(min([pix1(2) pix2(2) pix3(2) pix4(2)])); %取得x方向的负轴超出的偏移量
for i=1-delta_y:height-delta_y
for j=1-delta_x:width-delta_x
pix=[i j 1]/rot; %用变换后图像的点的坐标去寻找原图像点的坐标，
%否则有些变换后的图像的像素点无法完全填充
float_Y=pix(1)-floor(pix(1));
float_X=pix(2)-floor(pix(2));
if pix(1)>=1 && pix(2)>=1 && pix(1) <= h && pix(2) <= w
pix_up_left=[floor(pix(1)) floor(pix(2))]; %四个相邻的点
pix_up_right=[floor(pix(1)) ceil(pix(2))];
pix_down_left=[ceil(pix(1)) floor(pix(2))];
pix_down_right=[ceil(pix(1)) ceil(pix(2))];
value_up_left=(1-float_X)*(1-float_Y); %计算临近四个点的权重
value_up_right=float_X*(1-float_Y);
value_down_left=(1-float_X)*float_Y;
value_down_right=float_X*float_Y;
imgn(i+delta_y,j+delta_x)=value_up_left*img(pix_up_left(1),pix_up_left(2))+ ...
value_up_right*img(pix_up_right(1),pix_up_right(2))+ ...
value_down_left*img(pix_down_left(1),pix_down_left(2))+ ...
value_down_right*img(pix_down_right(1),pix_down_right(2));
end
end
end
subplot(122)
imshow(uint8(imgn))
Write Image: imwrite()
Lecture 09 Image Processing 图像处理 II：
问题的提出：
rice.png 中有多少颗米，以及每粒米的大小？
Image thresholding （图像阈值化）
A gray-level image can be turned into a binary image by using a threshold
graythresh() and im2bw()
graythresh() 用来计算最优的threshold 值
im2bw() Convert image to binary image by thresholding. b and w in im2bw() is b: black 0 w: white 1
%% 对一个灰度图进行 threshold
I = imread('rice.png');
subplot(131);
imshow(I);
% 首先找出一个最优的threshold level 值,这里使用 graythresh()
threshold_level = graythresh(I);%LEVEL is a normalized intensity value that lies in the range [0, 1]
subplot(132);
imhist(I); % 直方图
threshold_val = threshold_level*255; % 阈值
ylim=get(gca,'Ylim'); % 获取当前图形的纵轴的范围
hold on;
plot([threshold_val,threshold_val],ylim);
subplot(133);
bw = im2bw(I,threshold_level);
imshow(bw);
产生的问题：
1， 在 Backgroud 上存在 一些白色的点 （这是因为图片上面比较亮，导致背景被误判前景 ）
2， 一些米 不见了 （这是因为图片下面比较暗，导致前景被误判背景 ）
如何解决这些问题呢？
1，要Estimate the background
2, 要Substract the background
Background image estimation （背景估计）
Estimation for the gray level of the background
%% background estimation
I = imread('rice.png');
subplot(131);
imshow(I);
BG = imopen(I,strel('disk',15));
subplot(132);
imshow(BG);
% BG in 3-d view
x = 1:8:size(BG,1);
y = 1:8:size(BG,2);
[X,Y] = meshgrid(x,y);
Z = BG(x,y);
ax3 = subplot(133);
surf(X,Y,Z); % 3d 图
colormap(ax3,'jet');
set(gca,'xlim',[1,256])% image is 256x256
set(gca,'ylim',[1,256])
set(gca,'zlim',[0,255])%gray value is [0,255]
xlabel('x'); ylabel('y');zlabel('z');
axis square;
可以看出背景的上面值高，下面值低。
Background image subtraction（背景剔除）
%% background 剔除
I = imread('rice.png');
subplot(131);
imshow(I);
BG = imopen(I,strel('disk',15));
subplot(132);
imshow(BG);
I2 = I - BG; % 这时 整体的图像就比较均匀了，就改善了 亮度 不统一的问题
subplot(133);
imshow(I2);
此时我们再次进行 threshold 就会好很多
%% thresholding on Background Removed Image
I = imread('rice.png');
subplot(221);
imshow(I);
subplot(222);
threshold_level = graythresh(I);
bw = im2bw(I,threshold_level);
imshow(bw);
subplot(223);
imshow(I);
subplot(224);
% 去除背景
BG = imopen(I,strel('disk',15));
I2 = I - BG;
threshold_level2 = graythresh(I2);
bw2 = im2bw(I2,threshold_level2);
imshow(bw2);
Component-connected labeling （连接元素的标定）
现在怎么算呢？可以用 Component-connected labeling 算法，
Matlab 为我们提供了一个函数 bwlabel()
%% bwlabel()
I=imread('rice.png');
BG=imopen(I, strel('disk', 15));
I2=imsubtract(I, BG);
threshold_level=graythresh(I2);
BW=im2bw(I2, threshold_level);
[labeled, numObjects]=bwlabel(BW, 8);% 参数 8是8连通（上下左右加四角）。也可以是4，4连通（上下左右）
% labeled
% 是标定的矩阵，它的大小和BW一致
% numObjects是总数
disp(numObjects)
现在我们知道了米粒的总数为 99 ，
%% What is the size of the largest grain? What is the mean size of the grains?
I = imread('rice.png');
BG = imopen(I,strel('disk',15));
I2 = imsubtract(I,BG);
level = graythresh(I2);
BW = im2bw(I2,level);
[labeled,numObjects] = bwlabel(BW,8);
[heigh,width] = size(labeled);
largest = 0;%用于存储最大米粒的像素个数
sum = 0;%用于计算平均数前的总数
for num = 1:1:numObjects %用于迭代米粒的个数标号
calc = 0; %初始化每粒米个数
for i = 1:1:heigh
for j = 1:1:width
if labeled(i,j) == num
calc = calc + 1;
end
end
end
if calc > largest
largest = calc;
end
sum = calc + sum;
end
disp(largest);
average = sum/numObjects;
disp(average);
Note： 参考代码都是网上随便找的，可能质量不是很好
label2rgb() 函数
可以用它来给 米粒 标上颜色
%% label2rgb()
I=imread('rice.png');
BG=imopen(I, strel('disk', 15));
I2=imsubtract(I, BG);
threshold_level=graythresh(I2);
BW=im2bw(I2, threshold_level);
[labeled, numObjects]=bwlabel(BW, 8);
RGB_label=label2rgb(labeled); imshow(RGB_label);
小练习：
如何变为color 图的红米呢？ 具体做法是在 r层 标红
I=imread('rice.png');
BG=imopen(I, strel('disk', 15));
I2=imsubtract(I, BG); level=graythresh(I2);
BW=im2bw(I2, level);
[labeled, numObjects]=bwlabel(BW, 8);
No=zeros(1,99);
for i=1:256;
for j=1:256;
for n=1:numObjects;
if labeled(i,j)==n;
No(n)=1+No(n);
else
end
end
end
end
max(No)
mean(No)
perfectnum=0;
for k=1:99;
if No(k)<80;
perfectnum=perfectnum+0;
else if No(k)>350;
perfectnum=perfectnum+2;
else perfectnum=perfectnum+1;
end
end
end
disp(perfectnum);
%以上为显示优化计算后的米的数量，下面显示柱状图与红米；
subplot(121);histogram(No);
RGB_label = label2rgb(BW,'flag','k');
subplot(122);imshow(RGB_label);
Note： 参考代码都是网上随便找的，可能质量不是很好
regionprops() 函数：
%% regionprops()
I=imread('rice.png');
BG=imopen(I, strel('disk', 15));
I2=imsubtract(I, BG);
threshold_level=graythresh(I2);
BW=im2bw(I2, threshold_level);
[labeled, numObjects]=bwlabel(BW, 8);
graindata = regionprops(labeled, 'basic');
disp(graindata);
OUTPUT：
99x1 struct array with fields:
Area
Centroid
BoundingBox
Interactive Selection: bwselect()
Lets you select objects using the mouse
%% bwselect
I=imread('rice.png');
BG=imopen(I, strel('disk', 15));
I2=imsubtract(I, BG);
threshold_level=graythresh(I2);
BW=im2bw(I2, threshold_level);
[labeled, numObjects]=bwlabel(BW, 8);
ObjI = bwselect(BW);
imshow(ObjI); % 用鼠标 左键选定后， 按 右键
Lecture 10 Differentiation and Integration（微分和积分） ：
微分的定义：
微分：differentitation 是 △ 的意思。 可以简单理解为对其在某一点求导。
Polynomial differentiation（多项式的微分）
polyval() 函数
%% polyval() 计算x 在 f多项式下的值 f = [9,-5,3,7]; % 多项式的系数 x = -2:1:5; % x 值 y = polyval(f,x);% y 值 。它按照 y = 9x^3 -5x^2+3x +7 计算y值 plot(x,y);
polyder() 函数
%% polyder() 求一个多项式的导数（微分） f = [9,-5,3,7]; % 多项式的系数 f_d = polyder(f);% f' 多项式的系数，即函数的导数 disp(f_d);
OUTPUT 27 -10 3
小练习： 给定 𝑓(𝑥) =9x^3 -5x^2+3x +7，求𝑓′(7)?
%% f(x) = 9x^3 -5x^2+3x +7 求 f'(7) f_d = polyder([9,-5,3,7]); % f'(x) 的系数 ret = polyval(f_d,[7,]); % 求f'(x) 这个多项式在 x = 7 的值 disp(ret);
OUTPUT 1256
小练习：
%% conv() 函数 Convolution and polynomial multiplication. 卷积与多项式乘法 （一维卷积）
% 所谓两个向量卷积，说白了就是多项式乘法
f1 = [5,-7,5,10]; % 第一个多项式的系数
f2 = [4,12,-3]; % 第二个多项式的系数
f = conv(f1,f2); % 一维卷积 f1*f2 之后 多项式的系数
x = -2:0.1:1; % x值
y = polyval(f,x); % y值
plot(x,y,'b--','linewidth',2);
hold on;
f_d = polyder(f);% f'(x) 的系数
y2 = polyval(f_d,x); % y' 值
plot(x,y2,'r-','linewidth',2);
legend('f(x)','f ''(x)'); % 输出单引号使用 转义字符无效，使用两个单引号即可
这里涉及到了卷积【扩展】，
参考文章： https://www.cnblogs.com/hyb221512/p/9276621.html 和 https://blog.csdn.net/leviopku/article/details/80327478
一维卷积 conv : 简单记就是两个多项式 相乘
二维卷积 conv2 或 多维卷积 convn : 分三种 full , same ,valid
以二维卷积展示三种模式：
full mode:
same mode:
valid mode:
%% conv2() 函数 二维卷积
X = ones([3,3]); % 3x3
K = ones([3,3]); % 3x3
disp('======full=======');
ret = conv2(X,K); % default is full mode
disp(ret);
disp('======same=======');
ret2 = conv2(X,K,'same'); % default is full mode
disp(ret2);
disp('======valid=======');
ret3 = conv2(X,K,'valid'); % default is full mode
disp(ret3);
OUTPUT
======full=======
1 2 3 2 1
2 4 6 4 2
3 6 9 6 3
2 4 6 4 2
1 2 3 2 1
======same=======
4 6 4
6 9 6
4 6 4
======valid=======
9
关于 卷积输出的大小 计算公式：
如滑动步长为1，图片大小为3x3，卷积核大小为3x3，卷积后图像大小：5x5
Polynomial Integration（多项式的积分）
polyint() 函数
%% polyint() 函数 f = [5,0,-2,0,1]; % f(x) 多项式的系数 f_int = polyint(f,3); % 常数项 constant 设置为3 disp(f_int); % 求 ∫ f(7)? ret = polyval(f_int,7); disp(ret);
Numerical differentiation（数值的微分）
diff() 函数：
%% diff() 函数
x = [1,2,5,2,1];
ret = diff(x); %calculates the differences between adjacent elements of a vector。
% 相邻的
disp(ret);
小练习：
%% Exercise: obtain the slope of a line between 2 points (1,5) and (2,7)
% 获取 斜率
x = [1,2];
y = [5,7];
slope = diff(y)/diff(x); % 斜率
disp(slope);
%%
x0 = pi/2;
h = 0.1;
x = [x0,x0+h];
y = sin(x);
res = diff(y)/diff(x);
disp(res);
How does ℎ affect accuracy?
%%
hs = [0.1,0.01,0.001,0.0001,0.00001,0.000001,0.0000001];
% f = @(x) x+1;
x0 = pi/2;
f = @(h) diff(sin([x0,x0+h]))/diff([x0,x0+h]);
for h = hs
ret = f(h) - cos(pi/2); % cos(pi/2) 是精确值
disp(ret);
end
OUTPUT （Error Value） -0.0500 -0.0050 -5.0000e-04 -5.0000e-05 -5.0000e-06 -5.0004e-07 -4.9960e-08
How to Find the 𝑓′ over An Interval [0,2𝜋]?
%% 前面是一个点的导数,这里是一个区间的导数
h = 0.5;
x = 0:h:2*pi;
y = sin(x);
ret = diff(y)./diff(x); % 多个值 除就要 加 . 了 ！
disp(ret);
hold on;
plot(x,y);
plot(x(1:size(x,2)-1),ret,'ro-'); % x 要少一个
legend('sin(x)','sin''(x)','location','best'); % 会自动选择一个最佳的位置
Second and Third Derivatives
%% 二次 和 三次 微分(导数)
h = 0.5;
x = 0:h:2*pi;
%y = x.^3;
y = sin(x);
y_1 = diff(y)./diff(x);
y_2 = diff(y_1)./diff(x(1:end-1));
y_3 = diff(y_2)./diff(x(1:end-2));
hold on;
plot(x,y,'r.-');
plot(x(1:end-1),y_1,'g*-');
plot(x(1:end-2),y_2,'b+-');
plot(x(1:end-3),y_3,'cx-');
legend('sin(x)','sin''(x)','sin''''(x)','sin''''''(x)','location','best');
Numerical and integration （数值的积分）
第一个是用矩形来做预测， 第二个是用梯形来做预测 。
Midpoint Rule Using sum()
%% ∫ 4x^3 dx [0,2] h = 0.05; x = 0:h:2; midpoint = ( x(1:end-1) + x(2:end) ) ./2; s = 4*midpoint.^3 * h; A = sum(s); disp(A); % 15.9950
Trapezoid Rule Using sum()
%% ∫ 4x^3 dx [0,2] trapezoid f = @(x) 4*x.^3; % x is a vector h = 0.05; x = 0:h:2; x1 = x(1:end-1); x2 = x(2:end); s = ( f(x1) + f(x2) ) /2 * h; A = sum(s); disp(A); % 16.0100
也可以使用trapz() 函数，
h = 0.05;
x = 0:h:2;
trapezoid = trapz(4*x.^3);% 它接收的参数是 每个点的y值
% 返回的是 sum ( f(x0)+f(x1) ) /2
A = trapezoid*h;
disp(A); % 16.0100
Second-order Rule: 1/3 Simpson’s
第三种积分方法，它一次积两个的面积。所以它要求 等分的数目 一定是偶数个！
%% 1/3 simpson's rule
h = 0.05;
x = 0:h:2; % 必须是 EVEN 段 ,也就是 EVEN + 1 个点。而且要取到 0 和 2
y = 4*x.^3;
A = h/3*(y(1) + 4*sum(y(2:2:end)) + 2*sum(y(3:2:end-2)) + y(end));
disp(A); % 16 精准！
Numerical Integration: integral() 函数
%% integral() 函数 y = @(x) sin(x); % 匿名函数 其实y 是个函数指针 res = integral(y,0,pi); disp(res); % 2.0000 y2 = @(x) x.^2; % 当然也可以对 多项式来积分 res2 = integral(y2,0,1); disp(res2); % 0.3333
Double and Triple Integrals （二重 三重积分）
%% 二重/ 三重积分 f1 = @(x,y) y.*sin(x) + x.*cos(x); res1 = integral2(f1,pi,2*pi,0,pi); disp(res1); % -3.5864 f2 = @(x,y,z) y.*sin(x) + z.*cos(y); res2 = integral3(f2,0,pi,0,1,-1,1); disp(res2); % 2.0000
Lecture11 Root Finding （方程式求根）
Problem Statement
Symbolic approach （符号逼近）（解析解）
Symbolic Root Finding Approach
%% 1 使用 syms 或者 sym() syms x; % 此时 x 就是一个符号了 ret = x + x + x; disp(ret); ret2 = ret/4; disp(ret2); %OUTPUT: %3*x %(3*x)/4
sym('x'); % 此时 x 就是一个符号了
ret = x + x + x;
disp(ret);
ret2 = ret/4;
disp(ret2);
%OUTPUT:
%3*x
%(3*x)/4
%% to define y = x^2 - 2x -8
x = sym('x');
y = x^2 - 2*x -8; % y 也自动是个 sym 了
disp(y);
Symbolic Root Finding: solve()
%% solve() y = x*sin(x) -x = 0
syms x; % x 是个符号
res = solve('x*sin(x) - x',x);
disp(res); % res 也是个符号
% 也可以如下
syms x;
y = x*sin(x) -x; % y 也是个符号
res2 = solve(y,x);
disp(res2);
%OUTPUT
% 0
% pi/2
%
% 0
% pi/2
Solving Multiple Equations
%% solve multiple equations syms x y; eq1 = x - 2*y - 5; eq2 = x + y - 6; res = solve(eq1,eq2,x,y); disp(res); % res 是个 结构体 disp(res.x); disp(res.y); %OUTPUT: % x: [1x1 sym] % y: [1x1 sym] % % x: [1x1 sym] % y: [1x1 sym] % %17/3 % %1/3
Solving Equations Expressed in Symbols
%% syms a b x; y = a*x^2 - b; res1 = solve(y); % 默认以 x 为 未知数 disp(res1); res2 = solve(y,a); % 以 a 为 未知数 disp(res2); res3 = solve(y,b); % 以 b 为 未知数 disp(res3); %OUTPUT % b^(1/2)/a^(1/2) % -b^(1/2)/a^(1/2) % %b/x^2 % %a*x^2
%% Exercise % 1 syms a b x y r; eq = (x-a)^2 + (y-b)^2 - r^2; res = solve(eq,x); disp(res); %OUTPUT % a + (b + r - y)^(1/2)*(r - b + y)^(1/2) % a - (b + r - y)^(1/2)*(r - b + y)^(1/2) % 2 syms a b c d; m = [a b;c d]; res2 = inv(m); % inv() Matrix inverse disp(res2); disp(m*res2); %OUTPUT %[ d/(a*d - b*c), -b/(a*d - b*c)] %[ -c/(a*d - b*c), a/(a*d - b*c)] % %[ (a*d)/(a*d - b*c) - (b*c)/(a*d - b*c), 0] % 没有化简 %[ 0, (a*d)/(a*d - b*c) - (b*c)/(a*d - b*c)]
Symbolic Differentiation and Integration （符号 微分和积分）
%% symbolic differentitation diff() syms x; y = 4*x^5; res = diff(y); disp(res); %OUTPUT: %20*x^4
%% Exercise
%1
syms x y;
f = exp(x^2) /(x^3 - x + 3);
res = diff(f,x);
disp(res);
%2
g = (x^2 + x*y - 1 )/(y^3 + x+ 3);
res2 = diff(g,x);
disp(res2);
%OUTPUT:
%(2*x*exp(x^2))/(x^3 - x + 3) - (exp(x^2)*(3*x^2 - 1))/(x^3 - x + 3)^2
%
%(2*x + y)/(y^3 + x + 3) - (x^2 + y*x - 1)/(y^3 + x + 3)^2
%% symbolic integration int() syms x; y = x^2*exp(x); z = int(y); disp(z); z = z-subs(z, x, 0); % subs 是 替代的意思，在z 中，x用0取代 disp(z); %OUTPUT %exp(x)*(x^2 - 2*x + 2) % %exp(x)*(x^2 - 2*x + 2) - 2
%% Exercise syms x; f = (x^2 - x+1) /(x+3); F = int(f); disp(F); res = subs(F,x,10) - subs(F,x,0); disp(res); %OUTPUT: %13*log(x + 3) - 4*x + x^2/2 %13*log(13) - 13*log(3) + 10
Numeric root solvers （数值根求解）
Symbolic vs. Numeric
fsolve() 函数
%% fsolve() 函数 f = @(x) 1.2*x + 0.3 + x*sin(x); ret = fsolve(f,0); % 初始根 （猜测一个数） disp(ret); %OUTPUT %Equation solved. % %fsolve completed because the vector of function values is near zero %as measured by the default value of the function tolerance, and %the problem appears regular as measured by the gradient. % %<stopping criteria details> % % -0.3500
%% Exercise f = @(x) [2*x(1) - x(2) - exp(-x(1));-x(1)+2*x(2)-exp(-x(2))]; res = fsolve(f,[-5,-5]); disp(res); %OUTPUT %Equation solved. % %fsolve completed because the vector of function values is near zero %as measured by the default value of the function tolerance, and %the problem appears regular as measured by the gradient. % %<stopping criteria details> % % 0.5671 0.5671
fzero() 函数
%% fzero() 函数 f = @(x) x.^2; res1 = fzero(f,0.2); % guess value is 0.2 disp(res1); res2 = fsolve(f,0.2); disp(res2); %OUTPUT: % NaN 因为没有 cross the x axis % % 0.0063
fsolve() fzero() 的options
%% fsolve() fzero() 的options
f = @(x) x.^2;
%tolerance 偏差
options = optimset('MaxIter',1e3,'TolFun',1e-10);
res = fsolve(f,0.1,options);
disp(res);
res2 = fzero(f,0.1,options);
disp(res2);
%OUTPUT:
%Equation solved.
%
%fsolve completed because the vector of function values is near zero
%as measured by the selected value of the function tolerance, and
%the problem appears regular as measured by the gradient.
%
%<stopping criteria details>
%
% 1.9532e-04
%
%Exiting fzero: aborting search for an interval containing a sign change
% because NaN or Inf function value encountered during search.
%(Function value at -1.37296e+154 is Inf.)
%Check function or try again with a different starting value.
% NaN
%%
f = @(x) x.^2 - 2; % cross the x axis
res3 = fzero(f,1.4); % guess value is 1.4
disp(res3);
res4 = fsolve(f,1.4);
disp(res4);
res3 = fzero(f,-1.4); % guess value is -1.4
disp(res3);
res4 = fsolve(f,-1.4);
disp(res4);
%OUTPUT:
%
% 1.4142
%
%
%Equation solved.
%
%fsolve completed because the vector of function values is near zero
%as measured by the default value of the function tolerance, and
%the problem appears regular as measured by the gradient.
%
%<stopping criteria details>
%
% 1.4142
%
% -1.4142
%
%
%Equation solved.
%
%fsolve completed because the vector of function values is near zero
%as measured by the default value of the function tolerance, and
%the problem appears regular as measured by the gradient.
%
%<stopping criteria details>
%
% -1.4142
Finding Roots of Polynomials: roots()
%% roots() 函数 only works for polynomials res = roots([1,-3.5,2.75,2.125,-3.875,1.25]); disp(res); %OUTPUT: % 2.0000 + 0.0000i % -1.0000 + 0.0000i % 1.0000 + 0.5000i % 1.0000 - 0.5000i % 0.5000 + 0.0000i
How Do These Solvers Find the Roots?（上面的fsolve 和 fzero 都是如何找到根的呢？）
第一个方法： Bisection Method（Bracketing）（二分法）
l : lower u:upper
第二个方法： Newton-Raphson Method (Open)
不收敛的情况：
The tangent lines of x^3 − 2x + 2 at 0 and 1 intersect the x-axis at 1 and 0 respectively, illustrating why Newton's method oscillates between these values for some starting points.
它一直在 1 和 0 之间震荡，这种时候应该采用 二分法 来解决。
Lecture 12 Linear equations and Linear system(线性方程 和 线性系统)：
Linear Equation（线性方程）
Why Matrix Form?
更方便计算，也更简洁
如何解线性方程组？
常用的两种方法：
1. Successive elimination (through factorization) 逐次消除（通过因子分解）
又可以分为 Gaussian Elimination 和 LU Factorization
a Gaussian Elimination
MATLAB 中使用 Gaussian Elimination – rref()函数
%% Gaussian Elimination rref() Reduced row echelon form.
A = [1,2,1;2,6,1;1,1,4];
b = [2;7;3];
res = rref([A,b]);
disp(res);
%OUTPUT:
% 1 0 0 -3
% 0 1 0 2
% 0 0 1 1
b LU Factorization lu() 函数
参考里面的LU分解部分： https://www.cnblogs.com/zach0812/p/13368110.html
%% LU Factorization
A = [1,1,1;2,3,5;4,6,8];
[L,U,P] = lu(A); % permutation(交换) matrix P so that P*A = L*U
% 交换矩阵 实现行交换 (因为主元位置出现0的情况需要进行 行交换 )
disp(L);
disp(U);
disp(P);
disp('====');
disp(P*A);
disp(L*U);
%OUTPUT：
% 1.0000 0 0
% 0.2500 1.0000 0
% 0.5000 0 1.0000
%
% 4.0000 6.0000 8.0000
% 0 -0.5000 -1.0000
% 0 0 1.0000
%
% 0 0 1
% 1 0 0
% 0 1 0
%
%====
% 4 6 8
% 1 1 1
% 2 3 5
%
% 4 6 8
% 1 1 1
% 2 3 5
Matrix Left Division: \ or mldivide()
Ax = b => x = A \ b;
Matlab 的左除 \ 原理是使用 Successive elimination 来计算，不是使用 Cramer’s method （克拉默法则）
2. Cramer’s method （克拉默法则）Ax = b , x =^-1*b
矩阵C是“代数余子式矩阵”（cofactor matrix）
C 的转置矩阵 C^T。即伴随矩阵（adjoint matrix）。 参考： https://www.cnblogs.com/zach0812/p/13368110.html
克莱姆法则 Cramer’s Rule x = A^(-1)b
现在我们已经知道了 A的求逆公式，那么 对于Ax = b ,那么
Example:
𝑥 + 2𝑦 + 𝑧 = 2
2𝑥 + 6𝑦 + 𝑧 = 7
𝑥 + 𝑦 + 4𝑧 = 3
%% Cramer's method inv() 函数 A = [1,2,1;2,6,1;1,1,4]; b = [2;7;3]; x = inv(A)*b; disp(x); %OUTPUT: % -3.0000 % 2.0000 % 1.0000
Cramer’s Method 的问题：
A的inverse 不一定存在，此时 A矩阵就是一个 singular （奇异）矩阵。
下面哪个矩阵更健康？
肯定是B，
%% cond() 函数，来计算 一个矩阵 k值 越小越好 A = [ 1 2 3; 2 4.0001 6; 9 8 7]; cond(A) B = [ 1 2 3; 2 5 6; 9 8 7]; cond(B) %OUTPUT %ans = % 4.3483e+05 %ans = % 45.5623
Linear System （线性系统）
Eigenvalues and Eigenvectors 特征值和特征向量
作用一： 简化Ab 的计算，将矩阵向量的乘法 转化为 向量的加法
作用二： 阐述系统A的放大，缩小作用
A 是一个系统，
，
当 λ > 1 的时候就是 系统 A就是 放大的作用
当 λ < 1 的时候就是 系统 A就是 缩小的作用
%% eig() 函数 A = [2,-12;1,-5]; [vec,diag] = eig(A); disp(vec); disp(diag); %OUTPUT: % 0.9701 0.9487 % 两个特征向量 [0.9701;0.2425] % 0.2425 0.3162 [0.9487;0.3162] % % -1 0 % 两个特征值 -1 % 0 -2 -2
Lecture 13 Statistics and Data_Analysis（统计与数据分析）
Main Statistical Methodologies
Descriptive Statistics（描述统计学）
Mean, Median, Mode, and Quartile
%% mean() 函数 求平均数 d = [1,3,5,5,5,5,7,9,9,9,10,13,14]; res = mean(d); disp(res); %OUTPUT: % 7.3077 %% median() 函数 求中位数 d = [1,3,5,5,5,5,7,9,9,9,10,13,14]; res = median(d); disp(res); %OUTPUT: % 7 %% mode() 函数 Most frequent values in array d = [1,3,5,5,5,5,7,9,9,9,10,13,14]; res = mode(d); disp(res); %OUTPUT: % 5
关于四分位数：
四分位数（英语：Quartile）是统计学中分位数的一种，即把所有数值由小到大排列并分成四等份，处于三个分割点位置的数值就是四分位数。
- 第一四分位数（
），又称较小四分位数，等于该样本中所有数值由小到大排列后第25%的数字。 - 第二四分位数（
），又称中位数，等于该样本中所有数值由小到大排列后第50%的数字。 - 第三四分位数（
），又称较大四分位数，等于该样本中所有数值由小到大排列后第75%的数字。
第三四分位数与第一四分位数的差距又称四分位距（InterQuartile Range, IQR）。
计算方法：
离散分布
对于离散分布，在选择四分位数值方面没有普遍的共识。
方法1
- 使用中位数将有序数据集分为两半。
- 如果原始排序数据集中的数据点数奇数，请不要在任何一半中包括中位数（排序表中的中心值）。
- 如果原始排序的数据集中有偶数个数据点，则将此数据集精确地分成两半。
- 下四分位值是数据下半部分的中位数。高四分位值是数据上半部分的中位数。
方法2
- 使用中位数将有序数据集分为两半。
- 如果原始排序数据集中的数据点数量奇数，请在两个半数中均包括中位数（排序表中的中心值）。
- 如果原始排序的数据集中有偶数个数据点，则将该数据集精确地分成两半。
- 下四分位值是数据下半部分的中位数。高四分位值是数据上半部分的中位数。
方法3 （Matlab 中使用）
- 如果数据点的数量为偶数，则方法3与以上任一方法相同。
- 如果有（4 n +1）个数据点，则较低的四分位数为第n个数据值的25％加上第（n +1）个数据值的75％；高四分位数是第（3 n +1）个数据点的75％加第（3 n +2）数据点的25％。
- 如果有（4 n +3）个数据点，则较低的四分位数是第（n +1）个数据值的75％加第（n +2）个数据值的25％；高四分位数是第（3 n +2）个数据点的25％加第（3 n +3）数据点的75％。
方法4
%% prctile() 函数 Percentiles of data set. 百分数
d = [6, 7, 15, 36, 39, 40, 41, 42, 43, 47, 49];
res = prctile(d,50); % 50
disp(res);
res = prctile(d,25); % 25
disp(res);
res = prctile(d,75); % 75
disp(res);
%OUTPUT
% 40
%
% 20.2500
%
% 42.7500
Range and Interquartile Range （范围和四分位间距）
%% range = max - min d = [6, 7, 15, 36, 39, 40, 41, 42, 43, 47, 49]; range = max(d) - min(d); disp(range); %% interquartile range = q3 - q1 q1 = prctile(d,25); % 25 q3 = prctile(d,75); % 75 interquartile_range = q3 - q1; disp(interquartile_range);
Variance and Standard Deviation （方差和标准差）
%% 方差 d = [6, 7, 15, 36, 39, 40, 41, 42, 43, 47, 49]; S = var(d); disp(S); %OUTPUT: % 251.9636 %% 标准差 d = [6, 7, 15, 36, 39, 40, 41, 42, 43, 47, 49]; s = std(d); disp(s); %OUTPUT: % 15.8734
Figures Are Always More Powerful
x = 1:14; freqy = [1 0 1 0 4 0 1 0 3 1 0 0 1 1]; subplot(1,3,1); bar(x,freqy); xlim([0 15]); subplot(1,3,2); area(x,freqy); xlim([0 15]); subplot(1,3,3); stem(x,freqy); xlim([0 15]);
Boxplot 箱形图
%% 箱形图 Boxplot marks = [80 81 81 84 88 92 92 94 96 97]; boxplot(marks) res = prctile(marks, [25 50 75]); disp(res); %OUTPUT: % 81 90 94
Skewness（偏度）
%% 偏度 skewness()
d = randn([100,3]); % 100x3 正态分布
idx= d(:,1) < 0;
d(idx,1) = 0; % 第一列 小于0 都赋值为0 众数和最小值都为0
d(d(:,3)>0,3) = 0; % 第三列大于0 的都赋值为0 众数和最大值都为0
boxplot(d, {'Right-skewed', 'Symmetric', 'Left-skewed'});
res = skewness(d);
disp(res);
%OUTPUT:
% 1.6558 0.2162 -1.9950
注： 偏度与平均值和中位数之间的关系不直接相关：偏度为负（左偏）的分布的均值可以大于或小于中位数，正偏度（右偏）也相同。
许多教科书教给了一个经验法则，指出均值在右偏斜下位于中间值的右侧，而在左偏斜下位于中间值的左侧。
该规则以令人惊讶的次数失败。
它可能在多峰分布或尾巴较长而另一尾巴较重的分布中失败。
例如，在美国家庭中成年居民的分布中，它是个右偏分布。但是，由于大多数情况都小于或等于众数，均值位于较重的左尾。结果，在右偏的情况下，平均值在中位数的右边的经验法则失败了。
Kurtosis（峰度）
%% kurtosis 峰度
d = randn([100,3])*3; % 100x3 正态分布
k1 = kurtosis(d(:,1));
k2 = kurtosis(d(:,2));
k3 = kurtosis(d(:,3));
disp(k1);
disp(k2);
disp(k3);
%OUTPUT:
% 2.2985
%
% 3.2461
%
% 2.9984
Inferential Statistics（推断统计学）
置信区间 和
假设检验
省略这部分！
Lecture 14 Curve Fitting and Interpolation(曲线拟合和插值)
Simple Linear Regression
这里的目的就是要找到两个 β0 和 β1 使得拟合的效果最好！
Linear Regression Formulation (线性回归公式)
需要对它 做微分！ （求导）
注意：这里只有 β0 和 β1 两个是未知的，
继续解：
MATLAB中具体可以使用polyfit() 来解，如下：
Polynomial Curve Fitting: polyfit()
%% polyfit() 函数 x =[-1.2 -0.5 0.3 0.9 1.8 2.6 3.0 3.5]; y =[-15.6 -8.5 2.2 4.5 6.6 8.2 8.9 10.0]; fit = polyfit(x,y,1); %第三个参数是1次的拟合 返回拟合好的两个参数 ，
disp(fit); % 将拟合好的结果 画出来
xfit = [x(1):0.1:x(end)];
yfit = fit(1)*xfit + fit(2);
plot(x,y,'ro',xfit,yfit); %(x,y) 画点 (xfit,yfit)画线
set(gca,'FontSize',14);
legend('data points','best-fit');
scatter() 和 corrcoef() 函数
%% scatter() and corrcoef()
x =[-1.2 -0.5 0.3 0.9 1.8 2.6 3.0 3.5];
y =[-15.6 -8.5 2.2 4.5 6.6 8.2 8.9 10.0];
scatter(x,y); % 画散点图
box on;
axis square;
r = corrcoef(x,y); % x 和 y 的相关系数
disp(r);
%OUTPUT:
% x y
% x 1.0000(x 和 x) 0.9202 (x 和 y)
% y 0.9202(y 和 x) 1.0000 (y 和 y)
%% corrcoef() 扩展
x = randn(6,1);
y = randn(6,1);
z = randn(6,1);
A = [x,y,z]; % x,y,z 分别是一个列向量
r = corrcoef(A); % x 和 y 和 z 的相关系数
disp(A);
disp(r);
%OUTPUT:
% -0.0450 0.6989 0.6261
% -0.1024 -0.0259 1.1436
% -0.3022 0.3281 -1.8217
% -0.8536 -0.7011 0.2969
% -0.7193 -2.1663 0.1585
% -0.9170 0.0498 0.0461
% x y z
% x 1.0000(xx) 0.5515(xy) 0.1103(xz)
% y 0.5515(yx) 1.0000(yy) -0.1167(yz)
% z 0.1103(zx) -0.1167(zy) 1.0000(zz)
不是越高越好，可能会过拟合，
Multiple Linear regression
上面是一个变量的回归，下面是多个变量
What If There Exists More Variables?
对两个变量的线性回归 直观的表示就要在三维的空间中，
%% multiple linear regression
load carsmall;
y = MPG;
x1 = Weight; x2 = Horsepower;
X = [ones(length(x1),1) x1 x2];
b = regress(y,X);
x1fit = min(x1):100:max(x1);
x2fit = min(x2):10:max(x2);
[X1FIT,X2FIT]=meshgrid(x1fit,x2fit);
YFIT=b(1)+b(2)*X1FIT+b(3)*X2FIT;
scatter3(x1,x2,y,'filled'); hold on;
mesh(X1FIT,X2FIT,YFIT); hold off;
xlabel('Weight');
ylabel('Horsepower');
zlabel('MPG'); view(50,10);
% [b,bint,r,rint,stats]=regress(y,X);
adding column of ones to make it suitable for matrix multiplication.
%% multiple linear regression load carsmall; y = MPG; x1 = Weight; x2 = Horsepower; X = [ones(length(x1),1) x1 x2 x1.^2 x2.^2 x1.*x2]; b = regress(y,X); x1fit = min(x1):100:max(x1); x2fit = min(x2):10:max(x2); [X1FIT,X2FIT]=meshgrid(x1fit,x2fit); YFIT=b(1)+b(2)*X1FIT+b(3)*X2FIT + b(4)*X1FIT.^2 + b(5)*X2FIT.^2 + b(6)*X1FIT.*X2FIT; scatter3(x1,x2,y,'filled'); hold on; mesh(X1FIT,X2FIT,YFIT); hold off; xlabel('Weight'); ylabel('Horsepower'); zlabel('MPG'); view(50,10); % [b,bint,r,rint,stats]=regress(y,X);
上面所涉及的都是线性回归，如果不是线性的呢？
可以使用cftool()
Interpolation（内插 ，插值）
回归的目的是 拟合整体的数据
插值的目的是 数据不够，需要估计出新的数据，即插值
In the mathematical field of numerical analysis, interpolation is a type of estimation, a method of constructing new data points within the range of a discrete set of known data points.
In engineering and science, one often has a number of data points, obtained by sampling or experimentation, which represent the values of a function for a limited number of values of the independent variable.
It is often required to interpolate, i.e., estimate the value of that function for an intermediate value of the independent variable.
Common Interpolation Approaches (常见的插值方法)
Linear Interpolation: interp1() 函数
%% interp1() 1-D interpolation X = 0:0.5:10; V = sin(X); hold on; plot(X,V,'ro-'); Xq = 0:0.3:3; % new Xq Vq = interp1(X,V,Xq); plot(Xq,Vq,'g*')
Spline Interpolation（样条插值）: spline()
%% spline() Cubic spline data interpolation. X = 0:0.5:10; V = sin(X); hold on; plot(X,V,'ro-'); Xq = 0:0.3:3; % new Xq Vq = spline(X,V,Xq); plot(Xq,Vq,'g*')
What Are Splines?
2D Interpolation: interp2() 函数
%% interp2() 2D Interpolation
x = -2:.5:2; y = -2:.5:3;
[X,Y] = meshgrid(x,y);
Z = X.*exp(-X.^2-Y.^2);
hold on;
surf(X,Y,Z);
plot3(X,Y,Z,'ro');
x_i = 0:0.3:2; % new x
y_i = 0:0.3:3; % new y
[X_i,Y_i] = meshgrid(x_i,y_i);
Z_i = interp2(X,Y,Z,X_i,Y_i); % new Z_i
plot3(X_i,Y_i,Z_i,'g*');
view(3); % 3-d view
2D Interpolation Using Spline
%% 2D Interpolation Using Spline
x = -2:.5:2; y = -2:.5:3;
[X,Y] = meshgrid(x,y);
Z = X.*exp(-X.^2-Y.^2);
x_i = 0:0.3:2; % new x
y_i = 0:0.3:3; % new y
[X_i,Y_i] = meshgrid(x_i,y_i);
Z_i = interp2(X,Y,Z,X_i,Y_i); % default is linear
Z_i_2 = interp2(X,Y,Z,X_i,Y_i,'spline'); % new Z_i
subplot(121);
hold on;
surf(X,Y,Z);
plot3(X,Y,Z,'ro');
plot3(X_i,Y_i,Z_i,'g*-');
view(3);
subplot(122);
hold on;
surf(X,Y,Z);
plot3(X,Y,Z,'ro');
plot3(X_i,Y_i,Z_i_2,'g*-');
view(3);
结束
最后编辑于 2021-01-18 12:40 历时16天 于家