1

我的数据集还有另一个问题。基本上，有一个具有相关特征的基因列表，包括位置麻木（第 3 和 4 列）和链方向（+ 或 -）。我正在尝试对位置进行计算，以使它们相对于每个基因的起始密码子类型（第二列），而不是整个基因组（就像现在一样）。问题是，计算只在 + STRAND 序列上执行， - STRAND 序列没有出现在输出中。下面是数据集、我的代码、输出以及我尝试过的示例。

这是数据集：

    GENE_ID TYPE    POS1    POS2    STRAND
PITG_00002  start_codon 10520   10522   -
PITG_00002  stop_codon  10097   10099   -
PITG_00002  exon    10474   10522   -
PITG_00002  CDS 10474   10522   -
PITG_00002  exon    10171   10433   -
PITG_00002  CDS 10171   10433   -
PITG_00002  exon    10097   10114   -
PITG_00002  CDS 10100   10114   -
PITG_00003  start_codon 38775   38777   +
PITG_00003  stop_codon  39069   39071   +
PITG_00003  exon    38775   39071   +
PITG_00003  CDS 38775   39068   +

这是代码：

import numpy
import pandas
import pandas as pd
import sys

sys.stdout = open("outtry2.txt", "w")
data = pd.read_csv('pinfestans-edited2.csv', sep='\t')
groups = data.groupby(['STRAND', 'GENE_ID'])

corrected = []

for (direction, gene_name), group in groups:
    ##print direction,gene_name
    if group.index[group.TYPE=='start_codon']:
        start_exon = group.index[group.TYPE=='exon'][0]
    if direction == '+':
        group['POSA'] = 1 + abs(group.POS1 - group.POS1[start_exon])
        group['POSB'] = 1 + abs(group.POS2 - group.POS1[start_exon])
    else:
        group['POSA'] = 1 - abs(group.POS2 - group.POS2[start_exon])
        group['POSB'] = 1 - abs(group.POS1 - group.POS2[start_exon])
    ##print group
    corrected.append(group)

以下是输出示例：

     + PITG_00003
    GENE_ID     TYPE         POS1   POS2   STRAND  POSA  POSB
8   PITG_00003  start_codon  38775  38777  +       1     3   
9   PITG_00003  stop_codon   39069  39071  +       295   297 
10  PITG_00003  exon         38775  39071  +       1     297 
11  PITG_00003  CDS          38775  39068  +       1     294

以前我得到一个数组值错误（制表符分隔的数据集 ValueError 具有多个元素的数组的真相是模棱两可的错误），但这已经得到了解决。所以接下来我尝试只做这部分：

import numpy
import pandas
import pandas as pd
import sys

##sys.stdout = open("outtry2.txt", "w")
data = pd.read_csv('pinfestans-edited2.csv', sep='\t')#,
              #converters={'STRAND': lambda s: s[0]})
groups = data.groupby(['STRAND', 'GENE_ID'])

corrected = []

for (direction, gene_name), group in groups:
    print direction,gene_name

并且输出打印出所有 GENE_ID 和它们的 STRAND 符号（+ 或 -），并且它对 + 和 - 序列都执行了此操作。所以在下面的某个地方，它没有在 STRAND 列中选择任何带有 - 的序列。

所以我尝试将其添加到原始代码中：

if direction == '+':
    group['POSA'] = 1 + abs(group.POS1 - group.POS1[start_exon])
    group['POSB'] = 1 + abs(group.POS2 - group.POS1[start_exon])
elif direction == '-':
    group['POSA'] = 1 - abs(group.POS2 - group.POS2[start_exon])
    group['POSB'] = 1 - abs(group.POS1 - group.POS2[start_exon])
else:
    break
print group
# put into the result array
corrected.append(group)

这是输出的最后，它打印了第一个 - 然后在结束前冻结了一段时间：

+
        GENE_ID     TYPE         POS1    POS2    STRAND  POSA  POSB
134991  PITG_23350  start_codon  161694  161696  +       516   518 
134992  PITG_23350  stop_codon   162135  162137  +       957   959 
134993  PITG_23350  exon         161179  162484  +       1     1306
134994  PITG_23350  CDS          161694  162134  +       516   956 
-
4

1 回答 1

2

这些行对我来说似乎很奇怪：

if group.index[group.TYPE=='start_codon']:
    start_exon = group.index[group.TYPE=='exon'][0]

第一个，我猜，只是试图检查该组是否具有起始密码子标记。但这没有意义，原因有两个。

(1) 如果只有一个 start_codon 条目并且它是第一个，那么条件实际上是错误的！

In [8]: group.TYPE == 'start_codon'
Out[8]: 
0     True
1    False
2    False
3    False
4    False
5    False
6    False
7    False
Name: TYPE

In [9]: group.index[group.TYPE == 'start_codon']
Out[9]: Int64Index([0], dtype=int64)

In [10]: bool(group.index[group.TYPE == 'start_codon'])
Out[10]: False

也许你想要any(group.TYPE == 'start_codon')，或(group.TYPE == 'start_codon').any()或sum(group.TYPE == 'start_codon') == 1或什么？但这也不对，因为

(2) 您的代码仅在start_exon设置时才有效。如果不是，那么它要么给出 aNameError要么回退到上次发生的任何值，并且你无法保证这将是一个合理的顺序。

如果我只是单独使用start_exon = group.index[group.TYPE=='exon'][0]，那么我会得到

In [28]: for c in corrected:
   ....:     print c
   ....:     
       GENE_ID         TYPE   POS1   POS2 STRAND  POSA  POSB
8   PITG_00003  start_codon  38775  38777      +     1     3
9   PITG_00003   stop_codon  39069  39071      +   295   297
10  PITG_00003         exon  38775  39071      +     1   297
11  PITG_00003          CDS  38775  39068      +     1   294
      GENE_ID         TYPE   POS1   POS2 STRAND  POSA  POSB
0  PITG_00002  start_codon  10520  10522      -     1    -1
1  PITG_00002   stop_codon  10097  10099      -  -422  -424
2  PITG_00002         exon  10474  10522      -     1   -47
3  PITG_00002          CDS  10474  10522      -     1   -47
4  PITG_00002         exon  10171  10433      -   -88  -350
5  PITG_00002          CDS  10171  10433      -   -88  -350
6  PITG_00002         exon  10097  10114      -  -407  -424
7  PITG_00002          CDS  10100  10114      -  -407  -421

我不知道这些值是否有意义，但它似乎并没有跳过任何东西。

于 2013-01-28T00:59:39.227 回答