Chapter C. An Introduction to Matlab

Copyright 2016 Jon Danielsson. Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except in compliance with the License. You may obtain a copy of the License at. http://www.apache.org/licenses/LICENSE-2.0. Unless required by applicable law or agreed to in writing, software distributed under the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the License for the specific language governing permissions and limitations under the License.

Listing C.1: Variables Last edited: 2011

x = 10;
y = 5;
z = x+y;
who
Your variables are:
ans  x    y    z


Listing C.2: Arrays Last edited: 2011

a = 1:5
a =
1     2     3     4     5
b = [1,2,3,4;5,6,7,8]
b =
1     2     3     4
5     6     7     8
size(b)
ans =
2     4


Listing C.3: Sequences Last edited: 2011

x = 0:2:10
x =
0     2     4     6     8    10
x = linspace(0, 10, 6)
x =
0     2     4     6     8    10
x = logspace(0,2,4)
x =
1.0000    4.6416   21.5443  100.0000


Listing C.4: Load files Last edited: 2011

sp=load('sp.dat');


Listing C.5: Download S&P-500 data in Matlab Last edited: 2011

price = hist_stock_data('01012000','31122000','^gspc');
y=diff(log(price.Close(end:-1:1)))
plot(y)


Listing C.6: Basic data analysis Last edited: 2011

length(y)
sum(y)
prod(y)
range(y)
mean(y)
median(y)
var(y)
std(y)
sqrt(var(y))
corrcoef(y)
skewness(y)
kurtosis(y)
quantile(y,0.01)
min(y)
max(y)
sort(y)
abs(y)
diff(y)


Listing C.7: Matrix calculations Last edited: 2011

A = [2,4;7,3]
B = [3,2;1,4]
A .* B
ans =
6     8
7    12
A * B
ans =
10    20
24    26
C = [1,2]
A * C
??? Error using ==> mtimes
Inner matrix dimensions must agree.
A * C'
ans =
10
13
diag(C)
ans =
1     0
0     2
eye(3)
ans =
1     0     0
0     1     0
0     0     1


Listing C.8: Useful matrix functions Last edited: 2011

M = [1 2;3 6]
rank(M)
inv(M)
det(M)
eig(M)
[V,D] = eig(M)
sqrtm(M)
trace(M)
null(M)
poly(M)
chol(M)


Listing C.9: Distribution functions Last edited: 2011

q = 0
p = 0.5
mu = 1
sigma = 2
normpdf(q,mu,sigma)
normcdf(q,mu,sigma)
norminv(p,mu,sigma)
normrnd(mu,sigma,10,2)
normfit(x)


Listing C.10: Common distributions Last edited: 2011

trnd(3)
chi2rnd(3)
unifrnd(1,2)
lognrnd(1,2)
exprnd(2)
gevrnd(10,1,2)
gprnd(10,1,2)
binornd(3,0.3)
poissrnd(0.4)


Listing C.11: Multivariate distributions Last edited: 2011

mvnrnd(mu,sigma)
mvtrnd(C,df)
copularnd('Gaussian',rho,N)
copularnd('t',rho,NU,N)
copularnd('Clayton',alpha,N)


Listing C.12: Testing for normality Last edited: 2011

alpha=0.01
[h,p,jbstat,critval] = jbtest(y,alpha)
h
p
jbstat = 20.0648
critval = 11.5875
[h,p,ksstat,critval]=kstest(y,[],alpha)
h = 1
p = 0
ksstat = 0.4814
critval = 0.1020
qqplot(y)
qqplot(x,y)


Listing C.13: Basic time series analysis Last edited: 2011

[ac,acstd] = sacf(y,20,1,0);


Listing C.14: A script M-file Last edited: 2011

a = randn(10);
b = a * 5 + 2;
plot(b)


Listing C.15: A function M-file Last edited: 2011

function  k = mykurtosis(x)
m4 = mean((x - mean(x)).^4);
k = m4/std(x)^4 - 3
end


Listing C.16: Calling function M-file Last edited: 2011

mykurtosis(y)
k = 1.3502


Listing C.17: Kurtosis function Last edited: 2011

function  k = mykurtosis1(x,a)
if nargin == 1
a = 3;
end
m4 = mean((x - mean(x)).^4);
k = m4/std(x)^4 - a


Listing C.18: A for loop Last edited: 2011

x = randn(1,5);
z = NaN(1,5);
for i=1:5
z(i)=x(i)+0.5;
end


Listing C.19: An if-else loop Last edited: 2011

a = 10;
if (rem(a,3))==0
disp('a is a multiple of 3')
else
disp('a is not a multiple of 3')
end


Listing C.20: A while loop Last edited: 2011

a = 1;
n = 1;
while a < 100
a = a*n;
n = n + 1;
end


Listing C.21: A switch-case loop Last edited: 2011

x = 30;
units='kg';
switch units
case{'kilograms','kg'}
y = x*1000;
case{'grams','g'}
y = x;
case{'ounce','oz'}
y = x*28.35;
case{'pounds','lb'}
y = x*453.6;
otherwise
disp(['Unknown Units:' units])
y = NaN;
end


Listing C.22: Loop avoidance Last edited: 2011

a = 1:1000
b = 1000 - a
ssum = 0;
for i = 1:1000
ssum = ssum + a(i)*b(i);
end
ssum = a*b'


Listing C.23: Maximum likelihood estimation Last edited: 2011

randn('seed',1);
N = 100;
x = randn(N,1)*2+3;
theta0 = [-2, 5];
global x;
[theta, likelihood_value] = fminunc(@likfunc, theta0)
theta = 2.8082    1.9366
likelihood_value = 116.0928


Listing C.24: Normal log likelihood Last edited: 2011

function loglik = likfunc(theta)
global x;
T = length(x);
mu = theta(1);
sigma2 = theta(2)^2;
loglik = 0.5*T*log(sigma2);
loglik = loglik +  0.5*(sum((x-mu).^2 / sigma2));