This quick start reference originates from the short matlab tutorial that I had compiled for my students at the University of Applied Sciences HSZT in 2004 for my lecture Data Warehouse and Data Mining. It was a quick introduction to matlab which effectively fulfilled its purpose because my students were in fact able to implement simple estimation and optimization models within a couple of days. I extended this short tutorial then with comparable R commands to make it useful to a broader audience.
Both matlab and R are for me excellent mathematic software. Both of them have all the tools and a library of mathematical functions to implement all kinds of statistical risk, estimation or optimization models I want. Both have interface functions for reading files and databases, as well as writing to them.
I generally prefer matlab to R because matlab has a more intuitive syntax and a more userfriendly development environment with excellent documentation. The syntax of matlab is closer to the generic language of math compared to R; R looks sometimes more like a programming language with all its object stuff and illegible object property names. Apart from these slight advantages that would appeal probably more to beginners, matlab has the Simulink which is a great tool for visual demonstrations of analytical models with nice inputoutput blocks.
matlab is available for students and academic staff at some universities, but it is quite expensive for commercial use (see www.Mathworks.com). R is opensource and free (see www.rproject.org), and it has a broad acceptance for valuation and risk modeling among many financial institutes. You can download and begin to use R within minutes. The recent versions of introduction and user reference documents written for R are also quite userfriendly.
You can download this quick start reference as a formatted pdf file from the download page together with the executable demonstration scripts for matlab and R. These demo scripts include all the commands included in this reference, and more. I would also recommend the very useful MATLAB/R reference by David Hiebeler also downloadable as pdf file.
Tunç Ali Kütükçüoglu, 14. March 2012
MATLAB

R

Useful help commands


% general help
>>help
% general help
>>help
% elementary math functions
>>help elfun
% how while works
>>help while
% search for keywords
>>lookfor optimization

# open general help document (online)
>help.start()
# get help for a specific function or word
>help(‘min’) # or
>?min
# search for keywords
>help.search(‘optimization’) # or
>??optimization

Current (working) directory


>>pwd
% You can easily change the working directory by using the directory toolbox in GUI

>> getwd()
# see menu “File/Change dir…” for changing the working directory

Comment symbol


%

#

Assignment


Assignment % Assignment with automatic output to console if there is no semicolon at the end of the statement
>> a = 5
a =
5
% Assignment without output to console
>> a = 5;

# Assignment
> A=5 # or
> A
> A
[1] 5
# R doesn’t display anything in the console automatically even if there is no semicolon at the end of the statement

Create vectors


>> v = [1 3 6 9]
v = 1 3 6 9 % Series: StartValue:Interval:EndValue
>> v = 1:5
v = 1 2 3 4 5 >> v = 1:2:9 v = 1 3 5 7 9 >> v = 10:1:5 v = 10 9 8 7 6 5 % vector with multiple (repeating) constant values
>> v = ones(1,4) * 2
v = 2 2 2 2 % horizontal and vertical vectors
>> w = [1 3]’
w = 1 3 % For matlab, a horizontal vector is a 1xN, a vertical vector is a Nx1 matrix. There is not a separate entity as “vector” in addition to “matrix”.

>v = c(1, 3, 6, 9)
> v [1] 1 3 6 9 # Series: seq(StartValue,EndValue,Interval)
>v = 1:5
> v [1] 1 2 3 4 5 > v = seq(1,9,2) > v [1] 1 3 5 7 9 > v = seq(10,5,1) > v [1] 10 9 8 7 6 5 # vector with multiple (repeating) constant values
> v = rep(3,4)
> v [1] 3 3 3 3 # horizontal and vertical vectors
# R doesn’t make a distinction between vertical and horizontal vectors. A “vector” is not a “matrix” for R.

Create matrices


% 2×3 matrix with given element values
% “;” as row delimiter
% “,” or blank as column delimiter
>> D = [1 2 3; 4 5 6]
D = 1 2 3 4 5 6 % Create 2×4 matrix with random elements (values uniformly distributed between 0 and 1)
>>A = rand(2,4)
A = 0.8147 0.1270 0.6324 0.2785 0.9058 0.9134 0.0975 0.5469 % 1×3 matrix with all 3s
>>B = ones(1,3) * 3;
% 3×2 matrix with all 0s
>>C = zeros(3,2);

# 2×3 matrix with given element values
# fill values rowwise
>D=matrix(c(1,2,3,4,5,6),
nrow=2,byrow=TRUE) >D [,1] [,2] [,3] [1,] 1 2 3 [2,] 4 5 6 # fill values columnwise
>D=matrix(c(1,4,2,5,3,6),nrow=2)
# Create 2×4 matrix with random elements (values uniformly distributed between 0 and 1)
>A = matrix(runif(2*4),2,4)
[,1] [,2] [,3] [,4] [1,] 0.2180 0.1598 0.7320 0.7478 [2,] 0.9897 0.5809 0.7781 0.7134 # 1×3 matrix with all 3s
>B = matrix(3,1,3)
# 3×2 matrix with all 0s
>C = matrix(0,3,2)

Show all predefined variables


>>who

>ls()

Size of matrix


% number of rows and columns
>> [row,col] = size(A);
% number of rows
>> row= size(A,1);
% number of columns
>> col= size(A,2);
% maximum dimension (max(#row, #col))
>>length(A);
>>length([1 2 3]); % total number of elements in a matrix
>>numel(A)

# number of rows
>row = nrow(A)
# number of columns
>col = ncol(A)
# maximum dimension (max(#row, #col))
>max(dim(A))
# total number of elements in a vector
>length(v)
# total number of elements in a matrix
>length(A)
# dimensions of a vector
> dim(c(1,2,3))
NULL 
Access elements of a vector or matrix


% element of a vector
>> a = v(2);
% access element of a matrix with row and column index
>> x = M(2,3);

# element of a vector
> v[2]
# access element of a matrix with row and column indexD
> x = D[2,3]

Transpose vectors and matrices


% make a horizontal vector (i.e. 1xN matrix) vertical, or vice versa
>> v = [1 2 3]’
v = 1 2 3 % transpose matrix
>> X = [1 2 3; 4 5 6]’
X = 1 4 2 5 3 6 
# Transpose operation converts a vector into a horizontal (1xN) matrix
> t(c(1,2,3))
[,1] [,2] [,3] [1,] 1 2 3 # transpose matrix
>D=matrix(c(1,4,2,5,3,6),nrow=2)
> D [,1] [,2] [,3] [1,] 1 2 3 [2,] 4 5 6 > t(D) [,1] [,2] [1,] 1 4 [2,] 2 5 [3,] 3 6 
Add column(s) or row(s) to a matrix


>> v = [5 10]’
v = 5 10 % horizontal (columnwise) concatenation
>> X = [1 2; 3 4];
>> X = [X, v] X = 1 2 5 3 4 10 % vertical (rowwise) concatenation
>>v = [5 10];
>> X = [1 2; 3 4]; >> X = [X; v] X = 1 2 3 4 5 10 
> v = matrix(c(5,10), nrow=2)
> v [,1] [1,] 5 [2,] 10 # horizontal (columnwise) concatenation
>X=matrix(c(1,2,3,4),nrow=2,
byrow=TRUE) >X = cbind(X,v) > X [,1] [,2] [,3] [1,] 1 2 5 [2,] 3 4 10 # vertical (rowwise) concatenation
> v = matrix(c(5,10), nrow=1)
>X=matrix(c(1,2,3,4),nrow=2,byrow=TRUE) >X=rbind(X,v) > X [,1] [,2] [1,] 1 2 [2,] 3 4 [3,] 5 10 
Matrix partitioning (submatrices)


>>X = [1 2 3 4; 5 6 7 8]
X = 1 2 3 4 5 6 7 8 % single row of a matrix
>> X(1,:)
ans = 1 2 3 4 % submatrix with specified rows and columns
>>S = X(1:2, 2:4)
S = 2 3 4 6 7 8 % submatrix with specified rows and columns
>>S = X([1 2], [2 4])
S = 2 4 6 8 
>X=matrix(c(1,2,3,4,5,6,7,8),
nrow=2,byrow=TRUE) > X [,1] [,2] [,3] [,4] [1,] 1 2 3 4 [2,] 5 6 7 8 # single row of a matrix
> X[1,] # returns a vector
[1] 1 2 3 4 >X[1,,drop=FALSE] [,1] [,2] [,3] [,4] [1,] 1 2 3 4 # single column of a matrix
>X[,2,drop=FALSE]
[,1] [1,] 2 [2,] 6 # submatrix with specified rows and columns
> S = X[c(1, 2), c(2, 4)]
> S [,1] [,2] [1,] 2 4 [2,] 6 8 
Insert a matrix into another matrix


>>X = [1 2 3 4; 5 6 7 8]
X = 1 2 3 4 5 6 7 8 >> M = X; % insert by replacing elements (element assignment)
>>Y = [10 11; 12 13];
>> M(:, [2 3]) = Y M = 1 10 11 4 5 12 13 8 % insert by (horizontal) extension
>>X = [X(:,1), Y, X(:, 2:4)]
X = 1 10 11 2 3 4 5 12 13 6 7 8 
>X=matrix(c(1,2,3,4,5,6,7,8),nrow=2,
byrow=TRUE) > X [,1] [,2] [,3] [,4] [1,] 1 2 3 4 [2,] 5 6 7 8 >M=X # insert by replacing elements (element assignment)
>Y=matrix(c(10,11,12,13),nrow=2,
byrow=TRUE) >M[,c(2,3)] = Y > M [,1] [,2] [,3] [,4] [1,] 1 10 11 4 [2,] 5 12 13 8 # insert by (horizontal) extension
>X = cbind(X[,1], Y, X[,2:4])
> X [,1] [,2] [,3] [,4] [,5] [,6] [1,] 1 10 11 2 3 4 [2,] 5 12 13 6 7 8 
Basic matrix operations


>>X = [2 4 9; 16 25 36];
% square root of all elements
>>C = sqrt(X)
C = 1.4142 2.0000 3.0000 4.0000 5.0000 6.0000 % square of all elements
>>C = C .^ 2;
% matrix matrix multiplication
>> v = [1 2 3]’
v = 1 2 3 >>Y = X * v Y = 37 174 % matrix scalar multiplication
>>0.1 * X
ans = 0.2000 0.4000 0.9000 1.6000 2.5000 3.6000 % matrix scalar division
>>X / 10
ans = 0.2000 0.4000 0.9000 1.6000 2.5000 3.6000 % elementwise multiplication >> X = [1 2 3 ; 4 5 6]
>>X .* X % elementwise division
>> X ./ X
ans = 1 1 1 1 1 1 % inverse matrix
>> X = [1 3; 2 4] ;
>> Z = inv(X) Z = 2.0000 1.5000 1.0000 0.5000 >> X * Z ans = 1 0 0 1 
>X = matrix(c(2,4,9,16,25,36), nrow=2,
byrow=TRUE) # square root of all elements
>C = sqrt(X)
> C [,1] [,2] [,3] [1,] 1.414214 2 3 [2,] 4.000000 5 6 # square of all elements
>C = C ^2
# matrix matrix multiplication %*%
>v = matrix(c(1,2,3), nrow=3)
> v [,1] [1,] 1 [2,] 2 [3,] 3 >Y = X %*% v > Y [,1] [1,] 37 [2,] 174 # matrix scalar multiplication
> 0.1 * X
[,1] [,2] [,3] [1,] 0.2 0.4 0.9 [2,] 1.6 2.5 3.6 # matrix scalar division
> X / 10
[,1] [,2] [,3] [1,] 0.2 0.4 0.9 [2,] 1.6 2.5 3.6 # elementwise multiplication
> X * X
# elementwise division
> X / X
[,1] [,2] [,3] [1,] 1 1 1 [2,] 1 1 1 # inverse matrix
>X = matrix(c(1,3,2,4), nrow=2,
byrow=TRUE) >Z = solve(X) > Z [,1] [,2] [1,] 2 1.5 [2,] 1 0.5 > X %*% Z [,1] [,2] [1,] 1 0 [2,] 0 1 
Solve a typical matrix equation: A x b = c –> b = ?


>>A = [1 2; 3 4];
>>c = [1 1]’ >>b = inv(A) * c b = 1.0000 1.0000 % Test the result
>>A * b
ans = 1.0000 1.0000 
>A=matrix(c(1,2,3,4),nrow=2,
byrow=TRUE) >c = matrix(c(1,1), nrow=2) >b = solve(A) %*% c > b [,1] [1,] 1 [2,] 1 # Test the result
> A %*% b
[,1] [1,] 1 [2,] 1 
Sorting vectors


>>v = [2 5 1 4];
% sort vector in ascending order
% w: sorted vector
% ind: element indices such that w = v(ind)
[w, ind] = sort(v,’ascend’)
w = 1 2 4 5 ind = 3 1 4 2 >> v(ind) ans = 1 2 4 5 
>v = c(2,5,1,4)
# sort vector in ascending order
# w: sorted vector
# ind: element indices such that w = v(ind)
>res=sort(v,index.return=TRUE)
>w = res$x >ind = res$ix > w [1] 1 2 4 5 > ind [1] 3 1 4 2 > v[ind] [1] 1 2 4 5 
Sorting matrices


>> X=[3 5 2 7;1 10 12 5; 3 2 8 10;
1 6 4 2] X =3 5 2 7 1 10 12 5 3 2 8 10 1 6 4 2 % sort rows after first column
>> sortrows(X, 1)
ans = 1 10 12 5 1 6 4 2 3 5 2 7 3 2 8 10 % sort rows after first column in descending order
>> sortrows(X, 1)
ans = 3 5 2 7 3 2 8 10 1 10 12 5 1 6 4 2 % sort rows after first and second columns
>> sortrows(X, [1 2])
ans = 1 6 4 2 1 10 12 5 3 2 8 10 3 5 2 7 % sort columns after the first row; sort rows of the transposed matrix
>> sortrows(X’,1)’
ans = 2 3 5 7 12 1 10 5 8 3 2 10 4 1 6 2 
> X = matrix(c(3,5,2,7,1,10,12,5,3,2,8,
10,1,6,4,2),nrow=4,byrow=TRUE) > X [,1] [,2] [,3] [,4] [1,] 3 5 2 7 [2,] 1 10 12 5 [3,] 3 2 8 10 [4,] 1 6 4 2 # sort rows after first column
> X[order(X[,1]),]
[,1] [,2] [,3] [,4] [1,] 1 10 12 5 [2,] 1 6 4 2 [3,] 3 5 2 7 [4,] 3 2 8 10 # sort rows after first column in descending order
> X[order(X[,1]),]
[,1] [,2] [,3] [,4] [1,] 3 5 2 7 [2,] 3 2 8 10 [3,] 1 10 12 5 [4,] 1 6 4 2 # sort rows after first and second columns
> X[order(X[,1],X[,2]),]
[,1] [,2] [,3] [,4] [1,] 1 6 4 2 [2,] 1 10 12 5 [3,] 3 2 8 10 [4,] 3 5 2 7 # sort columns after the first row
> X[,order(X[1,])]
[,1] [,2] [,3] [,4] [1,] 2 3 5 7 [2,] 12 1 10 5 [3,] 8 3 2 10 [4,] 4 1 6 2 
Aggregation functions like sum, mean, max, min, stdev, variance


>> X = [1 4 2; 3 6 2]
X = 1 4 2 3 6 2 % sum of each column (columnwise sum)
>>sum(X)
ans = 4 10 4 % Generally, func(X, dim):
% dim = 1 ‘ columnwise operation
% dim = 2 ‘ rowwise operation
% default dim (dimension) is 1 (columnwise) if no explicit dimension is given
% sum of each row
>>sum(X, 2)
ans = 7 11 % sum of all elements of a matrix
>>sum(sum(X));
% mean value of each column
>> mean(X)
ans = 2 5 2 % look for other functions like min(), max(), std(), var(), cov(), median() in matlab help

>X = matrix(c(1,4,2,3,6,2),
nrow=2,byrow=TRUE) > X [,1] [,2] [,3] [1,] 1 4 2 [2,] 3 6 2 # sum of each column (columnwise sum)
>colSums(X)
[1] 4 10 4 # sum of each row
>rowSums(X)
[1] 7 11 # Note that the functions colSums() and rowSums() returns vectors; not matrices. # sum of all elements of a matrix >sum(A)
# mean value of each column
>colMeans(X)
[1] 2 5 2 # mean value of each row
>rowMeans(X)
[1] 2.333333 3.666667 # see other functions like apply(X,1,sd), apply(X,1,var), cov() and median() in R

Text output to console


% display simple text only
>> disp(‘this is a sentence’);
this is a sentence % text concatenate
>> str = [‘motor ‘, ‘car’, ‘s’]
str = motor cars >> v = [1 2 3 4]; >> disp([‘v = ‘, num2str(v)]); v = 1 2 3 4 % text output with an integer parameter, “\n” for explicit line feed
>> fprintf(‘a is equal to %d\n’, 5)
a is equal to 5 % text output with a floating number
>>fprintf(‘b is equal to %f\n’,
sqrt(2)) b is equal to 1.414214 % number formatting: Show two decimals after fractional point
>> fprintf(‘b is equal to %.2f\n’,
sqrt(2)) b is equal to 1.41 >> fprintf(‘Results: a = %d, b = %.1f, c = %.1f\n’, 4, 5.18, 6.12); Results: a = 4, b = 5.2, c = 6.1 
# display simple text only
>print(‘this is a sentence’);
[1] “this is a sentence” # text concatenate
>str=paste(‘motor ‘,’car’,’s’,sep=”);
> str [1] “motor cars” >v = c(1,2,3,4) > print(paste(‘v =’, paste(as.character(v),collapse=’ ‘))) [1] “v = 1 2 3 4” # text output with an integer parameter
>str=sprintf(‘a is equal to %d’, 5)
> str [1] “a is equal to 5” # text output with a floating number
> str = sprintf(‘b is equal to %f’,
sqrt(2)) > str [1] “b is equal to 1.414214” # number formatting: Show two decimals after fractional point
> sprintf(‘b is equal to %.2f’,
sqrt(2)) [1] “b is equal to 1.41” > sprintf(‘Results: a = %d, b = %.1f, c = %.2f’, 4, 5.18, 6.12); [1] Results: a = 4, b = 5.2, c = 6.12 
Programming language constructs


% while loop
>> v = [8 5 3 7 2 8 1 2];
% start from left, find the first number smaller than or equal to 2
>> i = 1;
>> while v(i) > 2 i = i+1; end >> fprintf(‘Results: i = %d,v(i) = %d\n’, i,v(i)) Results: i = 5, v(i) = 2 % alternative method
>>a = find(v>i = a(1)
% look for find() in matlab help: It is a useful and versatile search function
% for loop
% generate identity matrix (all diagonal elements 1, others 0)
>> I = zeros(3,3);
>> for i=1:3 I(i, i) = 1; end >> I I = 1 0 0 0 1 0 0 0 1 % alternative method
>> I = eye(3);
% if statements
>> x = 5;
% check if x is in range (2,10)
>> if (x > 2) && (x < 10)
disp(‘x is in range’) else disp(‘x is not in range’) end x is in range % see keywords like switch, break and continue for other programming constructs

# while loop
> v = c(8,5,3,7,2,8,1,2)
# start from left, find the first number smaller than or equal to 2
> i=1
> while (v[i]>2){ i=i+1 } >str=sprintf(‘Results: i = %d, v[i] = %d’,i,v[i]) > print(str) [1] “Results: i = 5, v[i] = 2” # alternative method
> a = which(vi = a[1]
# look for which() in R help; it is a useful search function
# for loop
# generate identity matrix (all diagonal elements 1, others 0)
> I = matrix(0, 3, 3)
> for (i in 1:3) { I[i,i]=1 } > I [,1] [,2] [,3] [1,] 1 0 0 [2,] 0 1 0 [3,] 0 0 1 # alternative method
> I = diag(3)
# if statements
> x = 5;
# check if x is in range (2,10)
> if^{1} {
print(‘x is in range’) } else { print(‘x is not in range’)} [1] “x is in range” # see keywords like switch, ifelse, repeat and next for other programming constructs

Writing scripts and functions


% script “test1.m” (text file) in current (working) directory
% start script
disp(‘This is a test script’); x=5 % end script % running the script from console
>>test1
This is a test script x = 5 % function with multiple arguments and a single return value: text file “difference.m” in current directory:
function c = difference(a,b)
c = a – b; end % calling the function from console
>> d = difference(10,2)
d = 8 % function with multiple return values; text file sumdif.m in current directory:
function [s,d] = sumdif(a,b)
s = a + b; d = a – b; end % calling the function from console
>> [s,d] = sumdif(5,3)
s = 8 d = 2 
# script “test1.r” (text file) in current (working) directory (menu “File>New script”)
# change the working directory from menu “File>Change dir…” if necessary
# start script
print(‘This is a test script’); x = 5 print(sprintf(‘x = %d’, x)) # end script # running the script from console
> source(‘test1.r’)
[1] “This is a test script” [1] “x = 5” # function with multiple arguments and a single return value: any text file, for example “test1.r” in current directory
difference = function(a,b) {
return(a – b) } # calling the function from console
> source(‘test1.r’)
> d = difference(10,2) > d [1] 8 # function with multiple return values; any text file, for example “test1.r” in current directory
sumdif = function(a,b) {
s = a + b d = a – b return(list(s, d)) } # calling the function from console
> source(‘test1.r’)
> results = sumdif(5,3) > s = results[[1]] > s [1] 8 > d = results[[2]] > d [1] 2 
Plotting graphs


% plot xsquare function
>>x = 1:100;
>> y = x .^ 2; >> plot(x, y); >> title(‘y = xsquare’) >> xlabel(‘x’) >> ylabel(‘y’) >> figure % opens new diagram % plot two curves on the same graph
% first curve is red, second is blue
>>x1 = 1:100;
>>y1 = x1 .^ 2; >>x2 = x1; >>y2 = (x2 – 5) .^ 2; >> plot(x1, y1,’r’, x2, y2,’b’) % histogram
% random number with normal distribution
>> x = randn(1,10000);
>> hist(x) % random number with uniform distribution
>> x = rand(1,10000);
>> hist(x) % z = (x + 0.5*y – 5) ^ 2
>> x = 0:0.1:5;
>> y = 0:0.2:10; >> [X,Y] = meshgrid(x,y); >> Z = (X + 0.5*Y 5) .^2; % surface plot
>> surf(X,Y,Z)
% contour lines (level curves)
>> contour(X,Y,Z)

# plot xsquare function
> x=1:100;
> y = x * x > plot(x,y,type=’l’,xlab=’x’, ylab=’y’,main=’y = xsquare’) > dev.new() # opens new diagram # plot two curves on the same graph
# first curve is red, second is blue
> x1 = 1:100
> y1 = x ^2 > y2 = (x5)^2 > plot(x,y1,type=”l”,col=”red”) > lines(x,y2,col=”blue”) # histogram
# random number with normal distribution
> x = rnorm(10000, sd=4, mean=5);
> hist(x) # random number with uniform distribution
> x = runif(10000);
> hist(x) # z = (x + 0.5*y – 5) ^ 2
> x = seq(0, 5, 0.1)
> y = seq(0, 10, 0.2) > f = function(x,y) return((x + 0.5*y 5)^2); > z = outer(x,y,f) # surface plot
> persp(x,y,z)
# contour lines (level curves)
> contour(x,y,z)

Constrained optimization
Solving optimization (max/min) problems subject to given constraints (boundary conditions)


% general syntax:
% x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub) – fun: singleoutput objective function whose return is to be minimized – x0: start value for vector x – inequality constraint: Ax ≤ b – equality constraint: Aeq x = beq – lb: lower bound for x – hb: higher bound for x % Problem:
Minimize z = sin(x + 0.5y) subject to constraints:
x + y ≥ 5 and x – y ≤ 2
% STEP 1: Create objective function objfun; text file objfun.m in current directory:
% Note: x –> v(1), y –> v(2)
function z = objfun(v)
z = sin(v(1) + 0.5*v(2)); end % STEP 2: Formulate boundary conditions
% inequlity condition:
% (x +y) ≤ 5
% x – y ≤ 2
>>A = [1 1; 1 1];
>>b = [5 2]; % there is no equality condition
>> Aeq = []; % empty matrix
>> beq = []; % there are no lower or upper boundaries for x or y (v(1) or v(2))
>>lb = [];
>>ub = []; % set a starting point for vector v
>>v0 = [0 0];
% STEP 3: Run optimization function
>>v = fmincon(@objfun,v0,A,b,
Aeq,beq,lb,ub); >> v v = 2.9302 3.5663 % Test the result
% Value of objective function at this point v
>> objfun(v)
z = 1.0000 % Inequality condition 1: check if x + y ≥ 5
>> sum(v)
ans = 6.4965 % Inequality condition 2: check if x – y ≤ 2
>> v(1)v(2)
ans = 0.6361 
# general syntax:
# x = constrOptim(x0, fun, grad, A, b, …) – fun: singleoutput objective function whose return is to be minimized – x0: start value for vector x – inequality constraint: Ax ≥ b – grad: gradient function for fun. It can be “null” if it is not known. # Problem:
Minimize z = sin(x + 0.5y) subject to constraints:
x + y ≥ 5 and x – y ≤ 2
# STEP 1: Create objective function objfun; in any text file, for example test1.r in current directory:
# Note: x –> v(1), y –> v(2)
objfun = function(v) {
z = sin(v[1] + 0.5*v[2]) return(z) } # declaring the function in console
> source(‘test1.r’)
# STEP 2: Formulate boundary conditions
# x + y ≥ 5
# x – y ≤ 2 –> x + y ≥ 2
> A = matrix(c(1,1,1,1),2,2,
byrow=TRUE) > b=matrix(c(5,2),2,1) >grad = NULL # set a starting point for vector v
# Note: initial value should satisfy the boundary conditions
> v0 = c(3,4)
# STEP 3: Run optimization function
>result = constrOptim(v0, objfun,
grad, A, b) > result $par [1] 2.740897 3.942860 $value [1] 1 # Test the result
# Value of objective function at this point v
> v = result[[1]]
> objfun(v) [1] 1 # Inequality condition 1: check if x + y ≥ 5
> sum(v)
[1] 6.683757 # Inequality condition 2: check if x – y ≤ 2 > v[1]v[2]
[1] 1.201963 
 x > 2) && (x < 10 [↩]