Monday, October 18, 2010

MATLAB

MATLAB, short for the MATrix LABoratory, is a strange program. Throughout my bachelor's and a bit during my master's I would hear people complaining about the program but I never had the opportunity to use it. Looking back on it now I wish I had!

One of the two courses I am taking this term is a modeling class; modeling how pollutants move/travel by dispersion/advection/diffusion/etc. I'm not sure how worthwhile the actual class will be because I don't forsee myself needing to model pollutants in the future but MATLAB has the opportunity to be a great tool.

The only problem is that I don't understand what I do wrong when it doesn't work! Take the following for example:

r=zeros(1,5);
i=1:600;
ufs=uf(i:end);
ufs=[uf(i:(i-1))];
Ru(i)=mean(uf.*ufs)


Do you understand what that is saying? Me either. That's what my professor gave me for the answer but it's not one that MATLAB is recognizing. I'm just trying to do an autocorrelation! And what I have above is not correct.

I'm sure MATLAB will be a very useful tool...once I figure it out. Until then it's a big pain.

4 comments:

  1. I don't follow that code at all! For a start, there is a definition of r which is not used after it is defined. Then, there is reference to a vector called uf which you haven't defined and so on.

    I haven't done autocorrelation before so I looked it up. It appears that there are various ways to define it. None of them seem particularly difficult to implement in MATLAB and I've done it for a couple of specific definitions.

    Let me know how to contact you and we can discuss which definition of auto-correlation you are using. Once I have this I'll send you the code.

    Cheers,
    Mike

    ReplyDelete
  2. I suppose I should have included my whole code instead of a little bit of it. The r (auto correlation)... I guess I'm not sure how to use after I've defined it. uf is really equal to "uf = u-mean(u);" and u is equal to a lot of things.

    I had originally been trying to use something like:

    autocorr(uf);
    [ACF,Lags,Bounds]=autocorr(uf,[],2)

    but my professor said that he would prefer that the auto correlation be done by hand.

    ReplyDelete
  3. Hello again. I don't have access to the autocorr function because that's part of the Econometric Toolbox. MATLAB has another built in function for auto correlation called xcov which is part of the Signal Processing toolbox (The Mathworks love their toolboxes!) and that's a toolbox I have. I've written a highly simplified auto-correlation function as follows

    %%%%BEGIN mikes_autocor.m
    function [ ac ] = mikes_autocor(x,maxlag,normalisation)
    %MIKES_AUTOCOR - calculates auto-correlation of a row vector x
    %results correspond to output given by MATLAB's xcov function
    %but without the repetition
    %This algorithm is inefficient but it does the job OK for small numbers of
    %lags
    %A better solution would be to use FFTs

    x = x-mean(x);
    xlen = length(x);

    if( nargin <2 )
    maxlag = xlen;
    end

    if(~isnumeric(maxlag))
    error('maxlag should be a number');
    end

    ac=zeros(1,maxlag);

    %The following loop calculates a 'raw' auto-correlation
    %If you don't give a normalisation then this is what will be output
    for lag=0:maxlag-1
    ac(lag+1) = (x(1+lag:xlen)*x(1:xlen-lag)');
    end

    if(nargin==3)
    if(strcmp(normalisation,'coeff'))
    normfactor = ac(1);
    ac = ac./normfactor;
    end

    if(strcmp(normalisation,'biased'))
    ac = ac./xlen;
    end

    if(strcmp(normalisation,'unbiased'))
    m = 0:maxlag-1;
    ac = ac./(xlen-m);
    end
    end

    %flip the results so that it matches the MATLAB xcov function
    %i.e. biggest lag returned first
    %remove the following you want zero lag first
    ac = fliplr(ac);

    end
    %%%%END mikes_autocor.m

    Now, there seems to be several definitions of auto correlation and auto-covariance floating about. I've taken MATLAB's xcov function as my basis which does an autocorrelation after subtracting the mean.

    One of the main differences between definitions seems to be the 'normalization' factor. Different disciplines seem to use their own one and then wonder why on earth you'd use anything else. xcov covers all bases and I've tried to replicate this in my function. Here's some sample output, comparing my function with xcov

    >> x=[1 2 3 4];
    >> xcov(x)
    ans =
    -2.2500 -1.5000 1.2500 5.0000 1.2500 -1.5000 -2.2500
    >> mikes_autocor(x)
    ans =
    -2.2500 -1.5000 1.2500 5.0000
    >> xcov(x,'biased')
    ans =
    -0.5625 -0.3750 0.3125 1.2500 0.3125 -0.3750 -0.5625
    >> mikes_autocor(x,4,'biased')
    ans =
    -0.5625 -0.3750 0.3125 1.2500
    >> xcov(x,'unbiased')
    ans =
    -2.2500 -0.7500 0.4167 1.2500 0.4167 -0.7500 -2.2500
    >> mikes_autocor(x,4,'unbiased')
    ans =
    -2.2500 -0.7500 0.4167 1.2500
    >> xcov(x,'coeff')
    ans =
    -0.4500 -0.3000 0.2500 1.0000 0.2500 -0.3000 -0.4500
    >> mikes_autocor(x,4,'coeff')
    ans =
    -0.4500 -0.3000 0.2500 1.0000

    The meaning of biased,unbiased and so on can be found at
    http://www.mathworks.com/help/toolbox/signal/xcorr.html

    My routine is extremely inefficient in that it uses a direct method of calculating auto-correlations. For large numbers of lags, the fastest way appears to use fast fourier transforms but I didn't fancy going there :)

    Let me know if you can't get the code working or if I've made a mistake.

    Hope this helps,
    Mike

    ReplyDelete
  4. Strangely enough, I feel like crying.

    Hope you figure this out.

    ReplyDelete

Related Posts Plugin for WordPress, Blogger...