Matlab: seqlogo with uniform plot column heights

In Matlab, I want to make a seqlogo plot o开发者_开发问答f an amino acid sequence profile. But instead of scaling the heights of the plot columns by entropy, I want all the columns to be the same height.

I'm in the process of modifying the code from the answers to this question, but I wonder if there is a parameter to seqlogo or some other function that I have missed that will make the column heights uniform.

Alternatively, is there a statistical transformation I can apply to the sequence profile to hack the desired output? (column heights uniform, height of each letter linearly proportion to its probability in the seqprofile)

Probably the easiest way around this problem is to directly modify the code for the Bioinformatics Toolbox function SEQLOGO (if possible). In R2010b, you can do:

edit seqlogo

And the code for the function will be shown in the editor. Next, find the following lines (lines 267-284) and either comment them out or remove them entirely:

S_before = log2(nSymbols);
freqM(freqM == 0) = 1; % log2(1) = 0

% The uncertainty after the input at each position
S_after = -sum(log2(freqM).*freqM, 1);

if corrError
    % The number of sequences correction factor
    e_corr = (nSymbols -1)/(2* log(2) * numSeq);
    R = S_before - (S_after + e_corr);
else
    R = S_before - S_after;
end

nPos = (endPos - startPos) + 1;
for i =1:nPos
    wtM(:, i) = wtM(:, i) * R(i);
end

Then put this line in their place:

wtM = bsxfun(@times,wtM,log2(nSymbols)./sum(wtM));

You will probably want to save the file under a new name, like seqlogo_norm.m, so you can still use the original unmodified SEQLOGO function. Now you can create sequence profile plots with all the columns normalized to the same height. For example:

S = {'LSGGQRQRVAIARALAL',...      %# Sample amino acid sequence
     'LSGGEKQRVAIARALMN',...
     'LSGGQIQRVLLARALAA',...
     'LSGGERRRLEIACVLAL',...
     'FSGGEKKKNELWQMLAL',...
     'LSGGERRRLEIACVLAL'};
seqlogo_norm(S,'alphabet','aa');  %# Use the modified SEQLOGO function

OLD ANSWER:

I'm not sure how to transform the sequence profile information to get the desired output from the Bioinformatics Toolbox function SEQLOGO, but I can show you how to modify the alternative seqlogo_new.m that I wrote for my answer to the related question you linked to. If you change the line that initializes bitValues from this:

bitValues = W{2};

to this:

bitValues = bsxfun(@rdivide,W{2},sum(W{2}));

Then you should get each column scaled to a height of 1. For example:

S = {'ATTATAGCAAACTA',...  %# Sample sequence
     'AACATGCCAAAGTA',...
     'ATCATGCAAAAGGA'};
seqlogo_new(S);            %# After applying the above modification

For now, my workaround is to generate a bunch of fake sequences that match the sequence profile, then feed those sequences to http://weblogo.berkeley.edu/logo.cgi . Here is the code to make the fake sequences:

function flatFakeSeqsFromPwm(pwm, letterOrder, nSeqsToGen, outFilename)
%translates a pwm into a bunch of fake seqs with the same probabilities
%for use with http://weblogo.berkeley.edu/

%pwm should be a 4xn or a 20xn position weight matrix. Each col must sum to 1
%letterOrder = e.g. 'ARNDCQEGHILKMFPSTWYV' for my data
%nSeqsToGen should be >= the # of pixels tall you plan to make your chart

[height windowWidth] = size(pwm);
assert(height == length(letterOrder));
assert(isequal(abs(1-sum(pwm)) < 1.0e-10, ones(1, windowWidth))); %assert all cols of pwm sum to 1.0

fd = fopen(outFilename, 'w');

for i = 0:nSeqsToGen-1
    for seqPos = 1:windowWidth
        acc = 0; %accumulator
        idx = 0;
        while i/nSeqsToGen >= acc
            idx = idx + 1;
            acc = acc + pwm(idx, seqPos);
        end
        fprintf(fd, '%s', letterOrder(idx));
    end
    fprintf(fd, '\n');
end

fclose(fd);
end