如何使用逻辑索引找到 DNA 序列子集的补码？答案

【问题标题】：How can I find the complement of a subset of a DNA sequence using a logical index?如何使用逻辑索引找到 DNA 序列子集的补码？
【发布时间】：2016-12-21 18:27:39
【问题描述】：

我有一个 DNA 序列，例如它的长度是 m*4n：

B = 'GATTAACTACACTTGAGGCT...';

我还有一个实数向量 X = {xi, i = 1..m*4n}，并使用mod(X,1) 将它们保持在 [0,1] 范围内。例如：

X = [0.223 0.33 0.71 0.44 0.91 0.32 0.11 ....... m*4n];

然后我需要通过应用函数将X 转换为二进制向量：

f(x)={0  ,0 < X(i,j) ≤ 0.5;  1 ,0.5 < X(i,j) ≤ 1;)

根据先前值的输出将类似于X = [0010100 ....]。如果X(i,j)==1，则补足B(i,j)，否则不变。在这种情况下，补码是匹配的碱基对（即 A->T、C->G、G->C 和 T->A）。

这是我迄今为止尝试过的代码，但没有成功：

%%maping X chaotic sequence from real numbers to binary sequence using threshold function
 X = v(:,3); 
 X(257)=[];
 disp (X);
 mode (X,1);
 for i=1
    for j=1:256
 if ((X(i,j)> 0) && (X(i,j)<= .5))
     X(i,j) = 0;
 elseif ((X(i,j)> .5) && (X(i,j)<= 1)) 
     X(i,j) = 1;
 end
    end
 end
 disp(X);

如何正确执行索引和补充？

【问题讨论】：

你需要更具体一些，我什至不知道 [图像处理] 是如何进入这个的，我不知道你想要什么。即便如此，也没有人会为您编写代码，因此您需要有一个您提出问题的工作/非工作示例。否则你的问题太宽泛了。
我加入了 Andras 的评论。请至少显示一个您正在尝试做的小数字示例。我完全不清楚你在问什么。
@AndrasDeak，先生，我确实尝试过自己写这篇文章，当我失败时，我来到这里请该领域的专家在这个特定部分帮助我。这是我整个程序的一小部分。感谢您的帮助
您应该添加您尝试过的内容，并描述它是如何失败的。正如我所说，我们可以帮助您调试代码，但我们不会为您编写代码。如果这个小部分太硬，我不确定你会在其余部分有更多的运气。你不同意吗？
我什至不确定您是如何将 9 个 4 字符的字符串编码为 3x3 矩阵的。我可以算出你的阈值函数（虽然很高兴知道Z(i,j)>1 会发生什么），但我希望你告诉我们 DNA 矩阵的补码是什么，而不是必须查找它。样本输入和期望的输出是必不可少的。

标签： matlab indexing dna-sequence

【解决方案1】：

给定一个存储为字符数组的示例碱基对序列：

B = 'GATTAACT';

还有一个与B长度相同的数值样本向量：

X = [0.223 0.33 0.71 0.44 0.91 0.32 0.11 1.6];

那么有一个相当简单的解决方案...

首先，您使用mod 函数意味着您只想使用X 中每个值的小数部分。你会这样做：

>> X = mod(X, 1)
X =
    0.2230    0.3300    0.7100    0.4400    0.9100    0.3200    0.1100    0.6000

接下来，您应该阅读documentation on vectorization。它将告诉您，在 MATLAB 中的许多操作都可以避免 for 循环。特别是，对您的向量 X 应用逻辑测试可以这样完成：

>> index = (X > 0.5)
index =
    0   0   1   0   1   0   0   1

而index 现在是logical index，长度与X 相同，每个大于0.5 的值都为1（即真）。您现在想要在B 中获取与这些索引对应的字符，将它们更改为它们的补码，然后将它们放回B。您可以使用 MATLAB 中的一个小技巧来做到这一点，即当使用 as 索引时，一个字符将转换为其 ASCII 数值：

>> compMap = '';  % Initialize to an empty string
>> compMap('ACGT') = 'TGCA'
compMap =
                                                                T G   C            A

注意字符'TGCA' 被放置在compMap 的索引65、67、71 和84 中（即'ACGT' 的ASCII 值）。其余都是空白。现在，您只需执行以下操作即可将索引碱基对替换为它们的补码：

>> B(index) = compMap(B(index))
B =
GAATTACA

综合起来，解决方案如下：

B = '...';     % Whatever your sequence is
X = [...];     % Whatever your values are
compMap = '';
compMap('ACGT') = 'TGCA';      % Build a complement map
index = (mod(X, 1) > 0.5);     % Get your logical index
B(index) = compMap(B(index));  % Replace with complements

【讨论】：

感谢先生您的回复..如果我想将其应用于整个矩阵（如果 B 是矩阵而不是 DNA 序列的子集）我应该改变什么@gnovice
@M.A.Fathy：如果X 是一个与B 完全相同大小的矩阵，那么上面的代码仍然可以工作，无需更改。
事实上，X 与 B 的大小不一样...它的大小是 256*1 ，，，我怎样才能使它的大小成为 256*4(256) ，有什么解决方案??我需要取X的每4位：例如（1001）来补充B的每个像素，例如（ACGT）..将其应用于整个矩阵B ...任何帮助请！ @gnovice