LESSON 16: Logarithmic scales and growth rates

FOCUS QUESTION: How can I use logarithmic scales to understand rates of growth?

Contents

EXAMPLE 1: Load the fetal size data

    load toRump.txt;  % Load the crown-to-rump data
    load toHeel.txt;  % Load the crown-to-heel data

EXAMPLE 2: Define variables to make code easier to read

   rWeeks = toRump(:, 1);          % Weeks in the toRump measurements
   hWeeks = toHeel(:, 1);          % Weeks in the toHeel measurements
   rMass = toRump(:, 3)/1000;      % Mass in kg of toRump measurements
   hMass= toHeel(:, 3)/1000;       % Mass in kg of toHeel measurements

EXAMPLE 3: Plot fetal mass as a function of gestational week

   figure
   hold on
   plot(rWeeks, rMass, 'ko')      % Plot the toRump mass
   plot(hWeeks, hMass, 'rs')      % Plot the toHeel mass
   xlabel('Week');
   ylabel('Mass (kg)');
   title('Fetal growth by gestational age');
   legend({'Crown-to-rump', 'Crown-to-heel'}, 'Location', 'NorthWest')
   hold off%%

EXAMPLE 4: Plot fetal mass as a function of gestational week using a semilog scale

   figure
   semilogy(rWeeks, rMass, 'ko', hWeeks, hMass, 'rs') % Plot the toRump mass
   xlabel('Week');
   ylabel('Mass (kg)');
   title('Fetal growth by gestational age');
   legend({'Crown-to-rump', 'Crown-to-heel'}, 'Location', 'SouthEast')

Questions Answers
What does semilogy(x, y) do? The semilogy function plots x versus the log10(y).
Why is the y-axis labeled in kg rather than the log of kg? The semilogyfunction, chooses major tick marks that correspond to powers of 10 in the original data. The ticks marks are labeled with the values in the original data. Hence, the reader uses original units.
Why not use two calls to semilogy with hold on? The semilogyfunction does not work with hold on. The plot reverts to linear scale.

EXAMPLE 5: Estimate fetal mass growth rates using the slope from the linear fit

   pRump = polyfit(rWeeks, log10(rMass), 1); % Fit weeks 8-20
   pHeel = polyfit(hWeeks, log10(hMass), 1); % Fit weeks 20-43
   corRump = corr(rWeeks, log10(rMass));     % Find correlation
   corHeel = corr(hWeeks, log10(hMass));

EXAMPLE 6: Output the estimates of fetal mass growth rates and R^2 for the linear fits

   fprintf(['Fetal mass growth rates in g/week:\n' ...
            '\t%g weeks 8-20, R^2 = %g\n' ...
            '\t%g weeks 20-43, R^2 = %g\n'], ...
             pRump(1), corRump*corRump, pHeel(1), corHeel*corHeel);

Fetal mass growth rates in g/week:
	0.209119 weeks 8-20, R^2 = 0.971853
	0.0491613 weeks 20-43, R^2 = 0.972093

Questions Answers
What is R2? The R2 is a measure of how much of the data variance is accounted for by the model (after the mean is taken out). A value of R2 close to one indicates a good model.
Is R2 always the square of the correlation? In general, no. However, for linear models (i.e., fitting a straight line), the value of R2 is the square of the correlation.

EXAMPLE 7: Read in the world population data and define variables

   data = csvread('WorldPopulation.csv');
   year = data(:, 1);
   population = data(:, 2);

EXAMPLE 8: Plot the world population as a function of time

Create a new cell in which you type and execute:

   figure
   plot(year, population)
   xlabel('Year');
   ylabel('Population (in billions)')
   title('World population (10,000 BC to 2050)');

Questions Answers
How can I interpret this graph? It is really hard to tell much except that in recent times the population grew very fast. The problem is that both the x and y data vary over too large a range to look at on a single graph. A logarithmic transformation should help here.

EXAMPLE 9: Plot year vs. world population on a semilog scale

   figure
   semilogy(year, population)
   xlabel('Year');
   ylabel('Population (in billions)')
   title('World population (10,000 BC to 2050)');

Questions Answers
How can I interpret this graph? This graph plots the year versus the logarithm of the population. If the population were growing exponentially, the graph would look linear. Unfortunately, this resulting graph appears to be growing "exponentially", meaning that the population growth rate is faster than exponential. (Think googel-plex rather than googel.) You will have to apply another log to find out just how bad things are.

EXAMPLE 10: Focus on the population data for 1500-2050

    yr = year(27:end);                      % Pick out years 1500-2050
    pop = population(27:end).*1e+9;         % Convert to units of people so log of log works
    popGRate = diff(pop)./diff(yr);         % Annual growth rate
    mYr = (yr(1:end-1) + yr(2:end))./ 2;    % The yr midpoints for plotting growth rate
    perCapitaGR = popGRate./pop(1:(end-1)); % Growth rate per person (birth rate)

Questions Answers
What does the per capita growth rate of a polynomial look like? The per capita growth rate of a polynomial (power of x) decays to zero.
What does the per capita growth rate of an exponential look like? The per capita growth rate of an exponential should be constant (flat).
How can I interpret the annual per capita growth rate of this graph (the green line)? From about 1650 to 1900 the per capita growth rate was flat and had a value of about 0.005, indicating the population was growing exponentially. After 1900 the per capita growth rate started to climb (population explosion). There was a temporary plateau between 1925 and 1945 at a growth rate of 0.01, before it took off again. The per capita growth rate is now falling but it is still substantially above rates from earlier times.

EXAMPLE 11: Plot the year vs the log of the log of the world population

   figure
   plot(yr, log10(log10(pop)))
   xlabel('Year');
   ylabel('Log_{10}(Log_{10}(Population))')     % Think about label here
   title('World population (1500-2050)')

Questions Answers
Why did this plot the log of the log of population? The previous graphs showed evidence that the population grew faster than exponential. If the rate of the exponential is growing exponentially, the graph of the log10(log10(pop) will appear linear (i.e., we have a googel-plex rather than merely a googel). This indeed appears to be the case.
How is plot(x, log10(log10(y)) different from loglog(x, y)? The plot(x, log10(log10(y)) draws a graph of x versus the log10(log10(y)) using an ordinary scale. The loglog(x,y) draws a graph of log10(x) versus log10(y) with both axes using logarithmic scales.

EXAMPLE 12: Plot log of the log of population versus its growth rate

Create a new cell in which you type and execute:

   figure
   ax = plotyy(yr, log10(log10(pop)), mYr, perCapitaGR);
   title('World population (1500-2050)')
   box off

   % Set the left axis labels and tick marks manually using the ax(1) handle
   xlabel(ax(1), 'Year');
   ylabel(ax(1), 'Population in billions (loglog scale)')
   set(ax(1), 'YLimMode', 'manual', 'YLim', [log10(8), log10(10)]);
   set(ax(1), 'YTickMode', 'manual', 'YTickLabelMode', 'manual');
   set(ax(1), 'YTick', [log10(8),  log10(9), log10(10)]);
   set(ax(1), 'YTickLabel', [0.1, 1, 10]);

   % Set the right axis labels and tick marks manually using the ax(2) handle
   xlabel(ax(2), '');                             % Turn off x-axis label
   ylabel(ax(2), 'Annual growth rate per capita') % Label right y-axis
   set(ax(2), 'XTickMode', 'manual', 'XTick', []); % Turn off x ticks

This lesson was written by Kay A. Robbins of the University of Texas at San Antonio and last modified on 31-Dec-2010. Please contact krobbins@cs.utsa.edu with comments or suggestions.


Published with MATLAB® 7.14