MatLab vs. Excel

The engineer who attempted teaching me MatLab used a short script to demonstrate its matrix manipulation prowess. MatLab (shortened from Matrix Laboratory) is built from the ground up to do matrix mathematics. I wanted to translate his script to VBA, to see, 1. If it’s possible, and 2. How Excel compares. This is the original MatLab script:

% fern.m        
%
% This is an exercise with MatLab Graphics and Flow Control.  The results
% should be a graphic drawn by an Iterated Function System.
%
% It takes about 2-3 minutes to complete the drawing
%
figure          % Establish a graphic window
pause(5.);      % An opportunity to re-size the Figure window
hold on         % Points accumulate in the window as an animation
x=[-5 0 5];     % 3 Initial points size the output axes
y=[0 1 10];
plot(x,y,’.g’,’MarkerSize’,4)
z = [x(2);y(2)];
%
pause(3.);     %  Pause before drawing
%
%  4 Matrices are key to the IFS Fractal description
%
fernM1 = [0 0
    0 0.16];
fernM2 = [0.2 -0.26
    0.23 0.22];
fernM3 = [-0.15 0.28
    0.26 0.24];
fernM4 = [0.85 0.04
    -0.04 0.85];
%
% Here is the fern
% We use tic/toc timing to measure execution time
%
tic
for kk = 1:20000
    dd = rand();
    if dd LTE 0.025
        z = fernM1 * z;
        plot(z(1),z(2),’.g’,’MarkerSize’,4)
        pause(0.01);        % Pause values less than 10ms makes no difference
    elseif (dd GT 0.025) & (dd LTE 0.125)
        z = fernM2 * z + [0 ; 1.6];
        plot(z(1),z(2),’.g’,’MarkerSize’,4)
    elseif (dd GT 0.125) & (dd LTE 0.225)
        z = fernM3 * z + [0 ; 0.44];
        plot(z(1),z(2),’.g’,’MarkerSize’,4)
    else
        z = fernM4 * z + [0 ; 1.6];
        plot(z(1),z(2),’.g’,’MarkerSize’,4)
    end
       
end
display(toc)
tic
for kk = 1:50000
    dd = rand();
    if dd LTE 0.025
        z = fernM1 * z;
        plot(z(1),z(2),’.g’,’MarkerSize’,4)
    elseif (dd GT 0.025) & (dd LTE 0.125)
        z = fernM2 * z + [0 ; 1.6];
        plot(z(1),z(2),’.g’,’MarkerSize’,4)
    elseif (dd GT 0.125) & (dd LTE 0.225)
        z = fernM3 * z + [0 ; 0.44];
        plot(z(1),z(2),’.g’,’MarkerSize’,4)
    else
        z = fernM4 * z + [0 ; 1.6];
        plot(z(1),z(2),’.g’,’MarkerSize’,4)
    end

end
display(toc)
display(‘The fern is done’)

 
Above, GT is “Greater Than”, LTE is “Less Than or Equals” and &amp is an ampersand. Here is my translation.

Sub Fern()
   Dim x(1 To 1, 1 To 3) As Double
   Dim y(1 To 1, 1 To 3) As Double
   Dim z()     As Variant
   Dim FernM1(1 To 2, 1 To 2) As Double
   Dim FernM2(1 To 2, 1 To 2) As Double
   Dim FernM3(1 To 2, 1 To 2) As Double
   Dim FernM4(1 To 2, 1 To 2) As Double
   Dim kk As Long, dd As Double
   Dim ZXArray(1 To 32000, 1 To 1) As Double
   Dim ZYArray(1 To 32000, 1 To 1) As Double
   Dim Fern As New Chart

   ReDim z(1 To 2, 1 To 1) As Variant
   x(1, 1) = -5: x(1, 2) = 0: x(1, 3) = 5
   y(1, 1) = 0: y(1, 2) = 1: y(1, 3) = 10
   z(1, 1) = x(1, 2)
   z(2, 1) = y(1, 2)

   FernM1(1, 1) = 0: FernM1(1, 2) = 0
   FernM1(2, 1) = 0: FernM1(2, 2) = 0.16

   FernM2(1, 1) = 0.2: FernM2(1, 2) = -0.26
   FernM2(2, 1) = 0.23: FernM2(2, 2) = 0.22

   FernM3(1, 1) = -0.15: FernM3(1, 2) = 0.28
   FernM3(2, 1) = 0.26: FernM3(2, 2) = 0.24

   FernM4(1, 1) = 0.85: FernM4(1, 2) = 0.04
   FernM4(2, 1) = -0.04: FernM4(2, 2) = 0.85

   For kk = 1 To 32000
      dd = Rnd()
      If dd LTE 0.025 Then
         z = Application.WorksheetFunction.MMult(FernM1, z)
      ElseIf dd GT 0.025 And dd LTE 0.125 Then
         z = Application.WorksheetFunction.MMult(FernM2, z)
         z(2, 1) = z(2, 1) + 1.6
      Elseif dd GT 0.125 And dd LTE 0.225 Then
         z = Application.WorksheetFunction.MMult(FernM3, z)
         z(2, 1) = z(2, 1) + 0.44      
      Else
         z = Application.WorksheetFunction.MMult(FernM4, z)
         z(2, 1) = z(2, 1) + 1.6
      End If
      ZXArray(kk, 1) = z(1, 1)
      ZYArray(kk, 1) = z(2, 1)
   Next kk
   Worksheets(“Sheet1”).Range(“A1:A32000”) = ZXArray
   Worksheets(“Sheet1”).Range(“B1:B32000”) = ZYArray

   Application.ScreenUpdating = False
   Set Fern = Charts.Add
   With Fern
      .ChartType = xlXYScatter
      .SetSourceData Source:=Sheets(“Sheet1”).Range(“A1:B32000”), PlotBy:=xlColumns
      .Location Where:=xlLocationAsNewSheet, Name:=”Fern”
      .Legend.Delete
      With .SeriesCollection(1)
         .MarkerBackgroundColorIndex = 10
         .MarkerForegroundColorIndex = 10
         .MarkerStyle = xlDiamond
         .Smooth = False
         .MarkerSize = 2
         .Shadow = False
      End With
      With .PlotArea
         .Border.LineStyle = xlNone
         With .Interior
            .ColorIndex = 2
            .PatternColorIndex = 1
            .Pattern = xlSolid
         End With
      End With
      With .Axes(xlCategory)
         .HasMajorGridlines = False
         .HasMinorGridlines = False
         .MajorTickMark = xlNone
         .MinorTickMark = xlNone
         .TickLabelPosition = xlNone
      End With
      With .Axes(xlValue)
         .HasMajorGridlines = False
         .HasMinorGridlines = False
         .MajorTickMark = xlNone
         .MinorTickMark = xlNone
         .TickLabelPosition = xlNone
         .MinimumScale = 0
         .MaximumScale = 10
         .MinorUnitIsAuto = True
         .MajorUnit = 2
         .Crosses = xlAutomatic
         .ReversePlotOrder = False
         .ScaleType = xlLinear
         .DisplayUnit = xlNone
         .Border.LineStyle = xlNone
      End With
   End With
   Application.ScreenUpdating = True
   
End Sub

 
Again, GT is “Greater Than” and LTE is “Less Than or Equals”. Make the substitutions after you paste in the code to a new workbook. We have to do this or the WordPress parser treats them as HTML tags and doesn’t present the complete information, trying as it does to helpfully shorten everyone’s code.

Some observations: You don’t dimension variables in MatLab; and case matters, with x being a different variable from X. MatLab can plot at least 70,000 points, while Excel is limited to 32,000. The Excel limit is supposedly a series limit, but I could never plot two series totaling over 32K points. MatLab doesn’t need a spreadsheet to plot, so if there’s a way in Excel to plot values without going through a spreadsheet, would someone please so describe. I don’t think you can. As expected, the Sheet1 ranges are the same dimensions as the ZX and ZY arrays, which are at the series limits, and allows copying the arrays to the spreadsheet.

Matrix Multiplication in VBA requires dimensioning the VB arrays as if they were cell arrays so that the worksheet function MMULT() can be used. Z() is a variant re-dimed to the dimensions of the matrix multiplication, meeting the specifications MMULT() requires in a spreadsheet–two rows by one column. I can’t make the code run on a Mac, where we’re stuck at VBA5. The Mac won’t accept the assignment of the array.

And finally, while it’s less than half the points, it’s much less than half the time to run. It’s almost instantaneous on my machine. That may mean the MatLab script is not optimum. I’ll pass on along any observations MatLab veterans may make.

Give it a try. You’ll see why it was chosen as a demo. It’s worth the trouble. A picture is worth one kilo-word, and you’ll have it in a blink.

…mrt

Posted in Uncategorized

19 thoughts on “MatLab vs. Excel

  1. So where is the extra time used in MATLAB’s implementation? Is it slower in general or, for example, is it faster at generating the points and slower at plotting them?

  2. Anyone doing matrix maths and manipulation should really look at APL. Ken Iverson’s language has been around since the 60s; it is clear, concise and compact.

    Cheers
    Andy

  3. For loops in Matlab are notoriously slow. You are supposed to take advantage of the various matrix and vector structures and avoid for loops as much as possible. So you would do something like rand([50000,1]) outside of a for loop to get dd to be a vector with 50000 numbers in it. Then you could use dd(i) to access each element in it within the for loop. I haven’t taken the time to figure out a way to get rid of everything else in the for loops, but taking the rand out of the loop may improve the speed a little.

  4. First off, that’s really cool. I used to try and do things like this comparing excel and vba, but since my vba abilities weren’t anywhere near as good as my matlab abilities there was no comparison.

    Couple observations from my experience w/ Matlab:
    1) Matlab handles arrays better than Excel in every possible regard – from creating them to doing any sort of analysis or calculation, it is easier and faster in Matlab. Case in point – scalar to multidimensional calculations can be done directly in Matlab whereas you need to loop through each element individually in Excel

    2) You can dimension arrays in Matlab, which speeds things up quite a bit w/ large data sets

    I don’t have Matlab on my computer any more so I can’t check it, but I’m guessing the reason there’s such a discrepancy between the two in terms of time is that the Matlab version is plotting each point as it’s calculated (matlab equivalent of leaving screenupdating on) vs. VBA does all the calculations then plots the result.

  5. I agree with Bjacobowski. I took out the plotting routines and just saved the data as the iteration proceeded. Then I did the plot at the end. The second loop (50,000 iterations) went from about 28 seconds on my laptop to about 0.08 seconds. I didn’t try the VBA version.

    If I do this with both loops and take out the pauses, the whole thing takes about 0.1 seconds.

  6. Michael,

    The reason that your MATLAB code is ‘slow’ is that every single iteration of the loop that computes the fern you are updating the plot. It is generally not a good idea to do this. I modified your code so that it simply computes the fern, and then plots it (see below). Running on even a fairly old tablet from IBM (X61) in R2010b the times I get for the 20000 iterations and then 50000 are

    >> fern
    ans =
    0.1405
    ans =
    0.3483
    The fern is done

    code – apologies if the HTML gets mangled …

    % fern.m
    %
    % This is an exercise with MatLab Graphics and Flow Control. The results
    % should be a graphic drawn by an Iterated Function System.
    %
    % It takes about 2-3 minutes to complete the drawing
    %
    figure(1) % Establish a graphic window – use the same on each time
    clf % Clear it and set the axis
    axis([-5 5 0 10]);
    hold on % Points accumulate in the window as an animation
    % Initial values
    z = [0;1];
    %
    % 4 Matrices are key to the IFS Fractal description
    %
    fernM1 = [0 0 ; 0 0.16];
    fernM2 = [0.2 -0.26 ; 0.23 0.22];
    fernM3 = [-0.15 0.28 ; 0.26 0.24];
    fernM4 = [0.85 0.04 ; -0.04 0.85];
    %
    % Here is the fern
    % We use tic/toc timing to measure execution time
    %

    % How many points are we making this time? Allocate the array
    N = 20000;
    allZ = zeros(2, N);
    tic
    for kk = 1:N
    dd = rand();
    if dd 0.025) && (dd 0.125) && (dd <= 0.225)
    z = fernM3 * z + [0 ; 0.44];
    else
    z = fernM4 * z + [0 ; 1.6];
    end
    allZ(:, kk) = z;
    end
    display(toc)
    % Plot what we have so far …
    plot(allZ(1,:), allZ(2,:), ‘.g’,’MarkerSize’, 4)
    pause(2); % See the sparse fern for 2 seconds
    % Now do more points
    N = 50000;
    allZ = zeros(2, N);
    tic
    for kk = 1:N
    dd = rand();
    if dd 0.025) && (dd 0.125) && (dd <= 0.225)
    z = fernM3 * z + [0 ; 0.44];
    else
    z = fernM4 * z + [0 ; 1.6];
    end
    allZ(:, kk) = z;
    end
    display(toc)
    % plot the rest (hold on so
    plot(allZ(1,:), allZ(2,:), ‘.g’,’MarkerSize’, 4)
    display(‘The fern is done’)

  7. The following comments have been received from Doug Jenkins over at Newton Excel Bach (see http://newtonexcelbach.wordpress.com/about/). Doug is a structural engineer and uses Excel, as he says, for functions without a dollar sign. He is a frequent commentator here, but at the moment is DDoE-challenged.

    ——-

    Expect some flack from the Matlab users (but on the other hand no “real” Matlab user would be seen dead at a site called Daily Dose of Excel, so maybe not :))

    There’s no doubt that for real matrix work Matlab will handle far bigger matrices, and deal with them far quicker, than the built in Excel functions, but if you are wanting to do that sort of thing in Excel there are tools available which will give you Matlab like functionality. For instance see:
    http://newtonexcelbach.wordpress.com/2010/05/27/compiled-alglib-matrix-functions-for-excel/

    I took the liberty of playing with your code to make the following changes:
    Read the number of graph points, and the basic “fern” parameters from the spreadsheet.
    Check if the chart sheet exists
    Alternative code using variants and doubles (with timer for each!)
    VBA matrix multiply instead of using worksheetfunction.

    Interesting what small changes to the parameters do to the fern plot. Excel 2010 will handle up to 1
    million data points by the way (although the fern looks nicer with about 30000!).

    I’d appreciate it if you would post the comments at DDoE for me.

    ——-

    I’m pleased to do so. Doug will be putting up his code over at his place. The fern is traveling across the Pacific.

    …mrt

  8. Jos –

    Thank you, but to be clear, it’s not my Matlab code, I’m just the translator. I don’t have the imagination to conceive either the elegance of the math or the beauty of the result. The author is a PhD mathematician who shared his time with such as me. I’ll pass along the comments and perhaps he’ll join us.

    …mrt

  9. I’ve heard from the author. He’s not the sort to blog, but he did give me this wiki pointer:

    http://en.wikipedia.org/wiki/Barnsley%27s_fern

    Michael Barnsley has several patents on fractal compression. At the above link are a good explanation of the math behind the fern, and two examples of “mutant” ferns from different coefficients.

    …mrt

  10. > so if there’s a way in Excel to plot values without going through a spreadsheet, would someone please so describe.

    Define names:

    Names.Add “x”, ZXArray
    Names.Add “y”, ZYArray

    Create the chart based on a dummy series eg: [a1:b2], then update the series values

    .XValues = “Sheet1!x”
    .Values = “Sheet1!y”

    Now you should be able to remove all worksheets and just leave the chart sheet.

    Btw, nice work! Lori

  11. Interesting post,
    The biggest performance hog in the Matlab code is, as Jos points out, the updating of the plot.

    With some help from the gurus over in the comp.soft-sys.matlab newsgroup, here you have an Matlab optimized version courtesy of Bruno Luong.

    Using this code the average of 50 trials gives me 0.0893 seconds to calc 50000 points (this is on a Intel Q9550, 4Gb ram, Vista, Matlab R2009a).
    Now that should surely knock the socks of any VBA routine :-P
    (still haven’t got completely rid of the for loop though)

    % Initial values
    z = [0;1];
    %
    % 4 Matrices are key to the IFS Fractal description
    %
    fernM1 = [0 0 ; 0 0.16];
    fernM2 = [0.2 -0.26 ; 0.23 0.22];
    fernM3 = [-0.15 0.28 ; 0.26 0.24];
    fernM4 = [0.85 0.04 ; -0.04 0.85];
    N = 50000;
    allZ = zeros(2, N);
    dd = rand(1,N);
    [na, where] = histc(dd, [0 0.025 0.125 0.225 1]);
    fern = { fernM1,fernM2,fernM3,fernM4 };
    a = cat(2, [0; 0], [0 ; 1.6], [0 ; 0.44], [0 ; 1.6]);
    tic
    for kk = 1:N
    w = where(kk);
    z = fern{w}*z + a(:,w);
    allZ(:,kk) = z;
    end
    display(toc)
    plot(allZ(1,:), allZ(2,:), ‘.g’,’MarkerSize’, 4)
    disp(‘The fern is done’);

  12. Peder –

    Thank you, and thanks to Bruno also. I forwarded your script to the author. There are things in that script I know weren’t in the course.

    Lori –

    Thank you for the compliment, and thank you for the tip. That’s a nice trick.

    …mrt

  13. Those are curly apostrophes in Peder’s script, compliments of WordPress. I can’t speak to Matlab, but I know Octave, the Gnu Matlab, doesn’t like them.

Leave a Reply

Your email address will not be published. Required fields are marked *