Using Recursive function to calculate all possible peptide combinations
13 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
Hello,
I'm trying to write a recursive function in MATLAB with input 'n' terms that calculates all possible peptides formed with the amino acids used being 'ACDEFGHIKLMNPQRSTVYW'. So if there are n characters in the sequence, there are 20^n combinations that could be produced from those characters. I'm trying to write a recursive function that would incorporate this. Any ideas?
Thank you
0 Commenti
Risposte (3)
Walter Roberson
il 22 Set 2017
There is a general computation pattern I refer to as the odometer pattern
possible_values = {'Z_Z', -3, 'V7', 'G', 'C', 119, 'Sr', 'Q', '2J'};
num_diff_vals = length(possible_values);
NumPlaces = 5; %or as needed
num_possibilities = num_diff_vals^NumPlaces;
results = cell(num_possibilities, NumPlaces);
odo = zeros(1, NumPlaces);
for K = 1 : num_possibilities
results(K,:) = possible_values( odo + 1 ); %note: adding 1 to each position
for P = NumPlaces : -1 : 1
odo(P) = odo(P) + 1;
if odo(P) == num_diff_vals
odo(P) = 0;
else
break;
end
end
end
The variable odo increases by 1 in the final digit each time. If that would take it past the maximum value, then it is reset to the first value and the next digit to the left is incremented. If that digit would overflow, zero it and increment the next digit to the left, and so on. At any given time, the odometer can be used as a vector index into the array of permitted values such as peptide letters.
This general structure has the same effect as if you had NumPlaces nested loops.
In the special case where all of the valid values are single characters, then you can use a character vector instead of a cell array to store the possible values, and a character array to store the results, with a minor change to how you initialize the results array.
0 Commenti
James Tursa
il 22 Set 2017
First, for even a moderate value of n, the number of possible combinations will be extremely large. Generating them as a list could take a very, very, long time and will likely use more memory than you have. Even if you could generate the list, it would take a long, long, time to process it in any algorithm you are envisioning.
Second, using a recursive function to generate this extremely long list would probably be the absolute slowest way to do it (unless you also mixed in some hard disk access which would make it even slower).
What size n are you thinking of using? What are you intending to do with this downstream in your code?
6 Commenti
James Tursa
il 22 Set 2017
Modificato: James Tursa
il 22 Set 2017
So, if you want to stick with the loops, here are some suggestions:
For the first column to start with, you don't need a loop at all. Just set the initial result to be the AA as a column vector. I.e.,
peptide = AA(:);
Then in the next loop, since you already have that first column, run the loop from i=2:n instead of i=1:n
Inside that loop, you will need to replicate the current peptide matrix prior to tacking on your AA(j) characters. You need to replicate it by the number of elements of AA (that is, each current row of peptide gets each AA element tacked on). That k loop then needs to run over the rows of the peptide matrix, appending your AA(j) characters. There is no need to call aminoAcid again in a recursive fashion. E.g., making these changes:
function [peptide] = aminoAcid(n)
%n is the amount of amino acids in the protein or peptide. This can also be
%referenced as the amount of characters. The output is the sequence for the
%protein or peptide.
AA = 'ACDEFGHIKLMNPQRSTVYW';
peptide = AA(:); % <-- set the first column directly, no need for a loop
m = numel(AA); % <-- remember the number of elements of AA
for i=2:n % <-- start the loop at 2
p = size(peptide,1); % <-- remember the current number of rows of peptide
peptide = repmat(peptide,m,1); % <-- replicate peptide by row for each element of AA
z = 1; % <-- start the generic counter
for j=1:m % <-- for each element of AA
for k=1:p % <-- for each row of peptide (prior to replicating it)
peptide(z,i) = AA(j); % <-- set the new character
z = z + 1; % <-- go to next row of peptide
end
end
end
Of course, much of this could be vectorized. But I was trying to retain as much of your looping logic as I could. The main things you were missing were replicating peptide at each of the outer loop iterations, and then noting that one of your inner loops should be using the number of rows of peptide instead of length(AA) for the indexing limit.
Jeremy Huard
il 12 Giu 2024
Your code could then look like the following:
AA = "ACDEFGHIKLMNPQRSTVYW";
% extract individual letters
letters = extract(AA,lettersPattern(1))
% generate all sequences of length 1 to 3
N = 3;
% preallocate a string array to store the results
Ncombinations = sum(numel(letters).^(1:N));
allcombinations = strings(Ncombinations,1);
% generate combinations
for i=1:N
repeatedletters = repmat({letters},1,i);
allcombinations_temp = combinations(repeatedletters{:});
idx = sum(numel(letters).^(1:(i-1)))+1 : sum(numel(letters).^(1:i));
allcombinations(idx) = join(allcombinations_temp{:,:},"",2);
end
allcombinations(17:24)
0 Commenti
Vedere anche
Categorie
Scopri di più su Nucleotide Sequence Analysis in Help Center e File Exchange
Prodotti
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!