Partially Finished Answer Page



The tutorial answer page is not quite finished yet. If you are really unhappy about this, then you can give me an incentive to finish by sending me email.

Matlab Tutorial (Answer Sheet)

Written by:

Steve McKelvey
Mathematics and Computer Science
Saint Olaf College

for the Envision It! Workshop, April 12, 1997

MATLAB and Matrices

Tutorial Task:


  1. Fire up MATLAB on your computer can create two 3x5 matrices and name them tut1 and tut2. The matrix you call tut2 should contain no zeros.

    The following commands should do the trick for you:

    tut1 = [1 2 3 4 5;6 7 8 9 10; 0 9 0 5 3]
    
    tut2 = [15 14 13 12 11; 10 9 8 7 6; 5 4 3 2 1]
    
  2. Use the size command to verify the size of the matrices you created above.

    The command size(tut1) and size(tut2) should both yield something like:

    ANS   =
     
        3.    5.
    

Shortcuts for Matrix Creation

Tutorial Task:


  1. Compare the result of the following two commands. Explain the difference in output.
    a1 = zeros(10)
    
    a2 = zeros(10);
    

    The first command will display the entire 10x10 matrix of zeros that it creates. The second will work silently. The difference is in the semicolon at the end of the command.

  2. What vector is created by the expression 1:2:10?

    The vector [1 3 5 7 9] is created.

  3. Create a vector of all even integers between 5 and 25.

    The following MATLAB command should do the trick:

    evens=[6:2:24]
    

Arithmetic Operations on Matrices

Scalar Operations

Tutorial Activity:

  1. Create matrices fdj and abc as indicated above and verify that abc was created successfully by entering the MATLAB command:
    abc
    
    which should show you the contents of the matrix abc.

    The matrix abc should consist of the entries:

        3.    6.    9.
       15.   12.    9.
       18.   15.   24.
    
  2. Determine the result of the following scalar arithmetic operations:
    fdj+3
    fdj-6
    fdj/2.
    
    in each case the indicated operation should be performed on each entry of the matrix fdj.

    The matrix fdj+3 is

    4. 5.  6.
    8. 7.  6.
    9. 8. 11.
    
    The matrix fdj-6 is
    -5. -4. -3.
    -1. -2. -3.
     0. -1.  2.
    
    The matrix fdj/2 is
    0.5  1.0  1.5
    2.5  2.0  1.5
    3.0  2.5  4.0
    

Matrix Addition and Subtraction

Tutorial Task

  1. Determine and explain the results of these MATLAB commands:
    abc = 1:10;
    def = 5:14;
    ghi = 3*abc + def
    

    The first creates the vector abc = [1 2 3 4 5 6 7 8 9 10]

    The second creates the vector def = [5 6 7 8 9 10 11 12 13 14]

    The third creates the vector ghi = [ 8. 12. 16. 20. 24. 28. 32. 36. 40. 44.]

  2. Determine and explain the results of these MATLAB commands:
    abc = [1 2 3 4;5 6 7 8];
    def = [4 3 2 1;0 -1 -3 3];
    abc + def
    

    The matrices abc and def are both 2x4 matrices, so they can be added together. The componentwise addition results in the matrix

    ANS   =
     
        5.    5.    5.    5.
        5.    5.    4.   11.
    
    

Componentwise Matrix Multiplication and Division

Tutorial Task

Create three 3x2 matrices named A, B and C and compute A.*B + C and (A+B)./C. Try to determine the results before asking MATLAB to perform the computation.

Here is an example:

A = <1 2;3 4;5 6>;
B = <2 3;4 5;6 7>;
C = <1 1;2 2;-1 -1>;
A.*B + C
ANS   =
 
    3.    7.
   14.   22.
   29.   41.

Componentwise Exponentiation

Traditional Matrix Multiplication

Tutorial Tasks

  1. Why does the last MATLAB command below fail?
    abc = [1 2 3 4];
    def = [2 5 4 3];
    abc * def
    

    The operation requested is a traditional matrix operation, which would require the vector def to be a column vector in this situation. (NOTE: It would also work if abc were a column vector while def remains a row vector. Try it and see what happens.)

  2. Why does the last MATLAB command below work?
    abc = [1 2 3 4];
    def = [2 5 4 3];
    abc * def'
    

    The apostrophe causes MATLAB to take the transpose of the vector def causing it to be tranformed into a column vector.

  3. This exercise demonstrates a quick way to total the entries in a vector. Suppose we have an array A (row vector) with, say, 123 entries in it. The sum of these entries can be easily computed using the MATLAB command:
    A*ones(1,123)'
    
    which multiplies the vector A with the column (after the transpose) vector made up of 123 ones. The result of this product is the total of the 123 entries in A.

    To explore this technique create a row vector with 10 entries, and use the idea above to direct MATLAB to compute the sum of the entries.

    a=[1 2 3 4 5 6 7 6 5 4];
    a*ones(1,10)'
     
    ANS   =
     
       43.
    
  4. Use the idea above, along with the short cut method for creating vectors whose entries have constant stepwidths, to create a single MATLAB command which sums up all the integers from 1 to 100. (HINT: The answer should be 5050.)

    [1:100]*ones(1,100)'
    
    ANS  =
    
         5050.
    
  5. What is the sum of all the even integers between 1 and 1001?

    [2:2:1001]*ones(1,500)'
     
    ANS   =
     
         250500.
    

Functions Applied to Arrays and Matrices

Tutorial Task

Using the window under the Apple Menu identify three more MATLAB functions which are usually applied to numbers but can also be applied to matrices within MATLAB.

There are dozens of these. I am sure you found some.

Two Dimensional Plots

Tutorial Tasks

  1. Create a graph of the function y=x/(1+x^2) for values of x between -5 and 5. First divide this interval into ten equal slices for graphing, and then do the same exercise dividing the interval into 100 equal slices. What effect do you see on the resolution of the plot?

    The commands for the first plot are:

    x = -5:1:5;
    y = x ./ (1 + x.^2);
    plot(x,y)
    

    The commands for the second plot are:

    x = -5:0.1:5;
    y = x ./ (1 + x.^2);
    plot(x,y)
    

    You should see an imporovement in the resolution with the increased number of grid points but you may also notice a significant increase in the time taken to produce the graph.

  2. Create a plot of a baseball diamond in as much detail as you would like. If you want to include a pitcher's mound you might want to take a look at the next task.

    The plot below should create the baseball diamond along with the foul ball lines extending into the outfield. Other possible enhancements could include small squares for the bases, detailed marking of the baselines, running coach stations, the batter's box and pitchers mound, as well as curves depicting the separation of the dirt and grass in the infield, etc.

    foulx = [ -150 0 150];
    fouly = [ 150 0 150];
    baselinex = [0 40 0 40 0];
    baseliney = [0 40 80 40 0];
    plot(foulx,fouly,baselinex,baseliney)
    
  3. What shape does the following set of MATLAB commands generate?
    theta = 0 : 0.05 : 2*pi;
    hold on
    axis('square')
    plot(cos(theta),sin(theta))
    

    This set of commands creates a circle of radius one centered at the origin. The ideas here are barrowed from the notion of polar coordinates in calculus.

    When you are finished using the image, type hold off. To see the effect of the axis('square') statement, enter the MATLAB command axis('normal') and redraw the graph.

    You will notice that the circle appears to be stretched. This is because the x and y graphs are on different scales. The axis('square') command forces the graph to use the same scales on both axes.

  4. As a final example, try entering these MATLAB commands. This graph is usually part of a calculus course dealing with polar coordinates.
    theta = 0 : 0.1 : 2*pi;
    r = sin(3*t);
    plot(r .* cos(theta), r .* sin(theta))
    
    Note the componentwise multiplication used in this problem, indicated by the dot-asterisk (.*) notation.

    You should see a three petaled rose-like figure appear in the plot.

    As an optional exercise (I'm a mathematician after all) guess what happens if the sin(3*t) term is changed to sin(4*t). Then test your guess by creating the graph.

    You will see an eight (surprise!!) petaled rose-like figure. You may have to increase the resolution (decrease the stepsize 0.1) to get a good picture of the figure.

Three Dimensional Graphs

Tutorial Tasks

  1. Using the help function of MATLAB, investigate the mesh command.

    You are on your own for this one.

  2. The function z=cos(2*(x+y)) produces surface plots of calm, peaceful waves. Graph this function for values of x and y between -10 and 10. Use stepsizes on the order of 0.2 or 0.3 when creating your x and yvectors.

    The MATLAB commands for this task are:

    x=[-10:0.2:10];
    y=[-10:0.2:10];
    [X Y] = meshgrid(x,y);
    Z = cos(2*(X+Y));
    surf(X,Y,Z)
    

Age-Class Population Model

Tutorial Task

  1. Why is the following MATLAB statement relevant to the above model? What value does it compute for us?
    ones(1,6)*(a^10*p)
    

    This statement sums up the populations in each age class, giving us the total population over all age groups.

  2. What MATLAB expression would give us the total number of individuals, regardless of age, in the critter population ten years hence?

    This is a trick question, the answer is the expression above.

  3. How would you quickly create a column vector showing the proportion of the total population in each of the age classes? This kind of information is very interesting in contexts such as the future social security system. Statements such as, "there will be only three workers for every retiree in the year 2020." come from this kind of analysis.

    The idea here is to divide the population in each age class by the total population. This quotient is a vector giving the proportion of the total population in each age class. The MATLAB command to do this, assuming p is the vector containing the current age class structure, is:

    proportion=p./(ones(1,6)*p)
    
    If you want to do this some years down the road, say 10 years, then you can replace each occurence of the vector p with (a^10*p).




    Answers only provided this far so far.






  4. Suppose a fertility drug is found which can increase the fecundity of two year olds, but leads to early deaths for older animals. In particular, suppose the drug increases the fecundity of two year old individuals to 2.3 births per year while the survival rates for individuals in age class 3 and 4 decrease from 0.8 to 0.6. An obvious question is whether this will increase the population of critters or not.

    To investigate this question we must update the entries in the A matrix to reflect the hypothetical introduction of the new drug. Only three entries inthe matrix A need be changed, so rather than recreate A from scratch we will make the three individual changes.

    The first change is that the birth rate for two year olds must be changed to 2.3. Since the birth rate for two year olds is found in the first row, third column of the matrix A, this change can be effected by the MATLAB command:

    A(1,3)=2.3;
    

    Make the remaining two changes in A, changes which reflected the lowered survival rates for age classes 3 and 4.

    Begin with the same initial population (all 10's) as before. Does the population seem healthier (whatever that means) under the drug regime? How does the drug effect the age distribution of the population?

  5. For this exercise we go back to the original Leslie setup, ignoring the changes you may have made in the previous exercise.

    If fecundity rates remain unchanged, what survival rates are needed for our critter population to remain stable? What proportion of a stable population falls into each age class?

It is fun to see how populations evolve over time. Toward this end it is possible to create a three dimensional graph showing the evolution of a Leslie population over time by creating a matrix whose columns are "snapshots" of the population at various times under study.

Suppose P is the initial population (all 10's in our example). We will use the matrix TimePop to be the matrix whose columns will contain the population at given points in time. Since the first population we are interested in is the initial population, we will set the first column of TimePop to the initial population using the MATLAB command

TimePop=P;
We will use a second matrix called PP to hold the most recent population computed. Thus, PP should start off holding the initial population
PP=P;
The task now facing us is to generate the population vector for the next time period, add this column vector to the matrix TimePop and then do the same again, moving ahead yet another time step. This can be achieved by executing the following MATLAB commands repeatedly:
PP = A * PP; TimePop = [TimePop PP];
The first command moves the population vector PP forward one time step. The second command adds a new column to the TimePop matrix. This new column contains the most recently computed population.

Using the up-arrow key, repeat the line above 15 to 20 times, creating a TimePop matrix with many population snapshots.

Tutorial Task

    The evolution of the critter population over time can be viewed by creating a three dimensional surface plot of the array TimePop. To create this surface plot enter the following MATLAB commands and note the plot which is generated. Pay special attention to the creation of labels and the colorbar.

    surf(TimePop)
    view(0,90)
    colormap(jet)
    colorbar
    xlabel('Year')
    ylabel('Age Class')
    

The colon (:) notation we discussed earlier is a convenient way to create large arrays where the entries are evenly spaced. A second use of the notation is to extract parts of existing matrices. In particular, a single colon is a reference to every row or every column in a matrix, depending on its position.

For example, the notation A(2,:) refers to every column of the second row of the matrix A, or, said another way, A(2,:) is the entire second row of A.

Similarly, B(:,5) refers to the entire fifth column of the matrix B.

Tutorial Task

  1. Continuing with the Leslie population model discussed above, what is the result of this MATLAB command?
    plot(TimePop(3,:))
    
  2. What MATLAB command would produce a plot of the population during the sixth time step?
  3. What MATLAB command would produce a plot showing the populations during both the third and sixth time steps?

Using calculators it is difficult to efficiently deal with Leslie models consisting of more than about ten age classes. Using MATLAB, however, it is relatively easy to deal with populations with dozens of age classes. Human populations are generally modeled with 101 age classes, allowing for ages between 0 and 100. It is also much easier to move far into the future, spotting long term trends.

Modifying the Leslie Model (Optional)

In this section we discuss modifications to the Leslie Population Model which describe a slightly different situation than that modeled by Leslie. This section introduces no new MATLAB skills and is present for those who want to undertake the challenge of creating a new model. I will describe the mathematics of the new model, but will leave it entirely up to you to implement the model in MATLAB. A knowledge of traditional matrix multplication, the non-componentwise variety, will be helpful.

In the Leslie model only two things can happen to individuals in an age class; they can die or they can move onto the next age class. In the model presented here we allow for a third option, remaining alive and in the same age class.

Given that it is possible to remain in the same class for several time steps, it is probably best to leave off the word age when describing these classes.

As far as I know, this version of the Leslie model has no well-known name, so we get to name the model for the purposes of this tutorial. Feeling singularly uncreative at the moment, I choose to call it Model X.

Model X requires three sets of parameters. Denote the fecundity of class i with the symbol f(i), the proportion of class i moving on to the next class at each iteration with the symbol s(i) and the proportion of each class which survives but remains behind with the symbol r(i).

Proceeding in a fashion reminiscent of the Leslie model development, if we create a model with classes numbered from zero to N, and let P(i) denote the number of individuals in class i, then we can derive the following equations which give the values of the various P(i) for the next time step in terms of the P(i) values for the current time step.

Considering, first, the values of P(i) in the next time step for classes numbered from one to N we see:

Next period's P(i+1) = s(i)*P(i) + r(i+1)*P(i+1)

This equation tells us that next year's class i+1 is made up of those individuals from this year's class i who move up a class as well as those individuals from this year's class i+1 who survive but do not advance in class. (a.k.a. flunking?) This equation holds for values of i from zero to N-1 inclusive.

A careful look at this situation reveals that we have not yet specified the population of class zero in next year's population. Recall that members of class zero are either true newborns or previous members of class zero who do not advance for some reason.

The mathematical expression which gives the number of class zero residents in next year's population is given by:

Next period's P(0)= (r(0)+f(0))*P(0) + f(1)*P(1) + f(2)*P(2) + ... + f(N)*P(N)

Pay special attention to the first term of this equation. It does not follow the pattern of the other terms.

Model X has a structure similar to the Leslie model. In fact, if P is a column vector whose N+1 entries are the populations of each class at some point in time, then there is an (N+1)x(N+1) matrix A so that the population next year can be determined by calculating the traditional matrix product

A*P

Tutorial Task

  1. For a population with six classes, use your imagination to determine values for s(i), r(i) and f(i) for each class. You can use the Leslie model example for inspiration.
  2. Using the values found above, design and create in MATLAB the 6x6 matrix A which can be used to move from one time period to the next in the model.
  3. Perform all the Tutorial Tasks from the section on Leslie models using, instead, Model X.


http://www.stolaf.edu/people/mckelvey/envision.dir/mar97.tutorial.html
Last Update: 11 MAR 1997

Disclaimer

Disclaimer