Matching Pursuit for one-dimensional signals: Matlab code

Matching Pursuit (MP) is a greedy algorithm to obtain a sparse representation for a signal in terms of elements of a dictionary. I’ve seen it used from time to time in neuroscience. For example, Smith and Lewicki (2006) use it as part of their demonstration that a sparse code for natural sounds matches the properties of cochlear receptors. It’s also frequently used in EEG/MEG analysis and recently was applied to LFPs in Ray and Maunsell (2011).

With one-dimensional signals such as EEG and sound, it’s natural to use a dictionary composed of a set of basis elements offset at every possible time point. In that case, the result is a sparse convolutional code.

Here’s an implementation in Matlab. It’s fairly efficient; it precomputes all convolutions outside the main loop and uses a 3-level search tree to find the current maximal correlation. I’ve also implemented a variant which forces all basis functions to be positive. It was borne out of my frustration in failing to figure out how to use the “anywave” block in the Matching pursuit toolkit (MPTK).

15 responses to “Matching Pursuit for one-dimensional signals: Matlab code”

  1. HELLO,I am learning your matching pursuit code ,I have some questions for you?Can you answer them for me?
    1.When build the gabors atoms,what’s the meaning of sigmas and rg?How to choose them?And how to build the atoms accroding my signal?
    2The output ‘ws’ ,can I use it for WVD?Can you give me some tips?

  2. Thank you for your reply,but I tried you suggestion and revised the code,the following is my code:
    and i got new r1 and r ,but there is stil big difference between them!Do you have details about the output ws and r,or example? Can you leave a contact way so that we can have a further discussion.
    rg = (-500:500)’;
    sigmas = exp(log(2):.3:log(200));
    gabors = exp(-.5*(rg.^2)*sigmas.^(-2)).*cos(rg*sigmas.^(-1)*2*pi*2);
    r1=0;
    [ws,r] = temporalMP(X,gabors,false,1000);%X is the singal to process,It ‘s size is 1000*1single
    gabors = bsxfun(@times,gabors,1./sqrt(sum(gabors.^2)));
    for i=1:16
    r1 = r1 + conv(ws(:,i),gabors(:,i),’same’);
    end

    • Now ,i have some ws and r ,if i want to do WVD time-frequnecy ,how to express ws and gabors’sTFR?

  3. Hello,Now I am use your temporalMP.m.I have you a question what’s the expression between the output ws and r? You wrote r = sum_i conv(ws(:,i),B(:,i),’same’),I tested it,but it is wrong.IF i have creat 16
    gabors ,my code is
    ;for i=1:16
    r1=r1+conv(ws(:,i),gabors(:,i),’same’);
    end
    I use it to test the relations between the r1 and r,but the outcome is aweful.Can you explain it for me?
    How to express the relations between ws and r?

    • It might be a scaling thing. Try using:

      gabors = bsxfun(@times,gabors,1./sqrt(sum(gabors.^2)));
      Before running your code. Also, if your dictionary elements are assymetric, you will also need to use:

      r1 = r1 + conv(ws(:,i),gabors(end:-1:1,i),’same’);

      During the reconstruction.

    • Hey Chen, I checked and the correct scaling is:

      gabors = bsxfun(@times,gabors,1./(sum(gabors.^2)));

      Sorry about that

  4. Hi, to understand the MP algorithm, I also want to use the “mpp” package from ftp://cs.nyu.edu/, when I use

    make -f make_linux

    on my ubuntu, several warnings come out , and indicates

    “mv main.o obj/ ”

    I did not know what should I do next, will plz help me? thanks.

  5. Hi, Patrick ,nice to find your matlab code of the matching pursuit algorithm, as a new comer in this field, I think your work will help me a lot, my problem is—

    when run the TryTemporalMP.m, the error just like “Unexpected MATLAB operator”

    Error: File: temporalMP.m Line: 172 Column: 11
    Expression or statement is incorrect–possibly unbalanced (, {, or [.
    e.g. [~,maxoffset] = max(Cmax(rg));

    Could you tell me how to solve this problem, by the way, my MATLAB is 7.5.0

    thanks.

    • Replace [~,maxoffset] = max(Cmax(rg)); with [dummy,maxoffset] = max(Cmax(rg)); ~ is only supported in later Matlab versions.

Leave a comment

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s