-2

我认为我的问题并不完全清楚。所以，是时候用另一种方法来解释我的追求了。

  • 我有一个包含大约 27,500 个数据点的数据文件。每个数据点在第一列中都有一个介于 0 和 46041637 之间的唯一整数，在第 2 - 6 列中有描述。
  • 我想知道整数如何分布在 46M 的可能性中。例如，有多少（以及哪些）数据点在 1-1000 范围内，有多少在 1001-2000 之间等等，直到达到 46M。
  • 在前面的示例中，将创建大约 46K (46M/1000) 个较小的数据集（或 bin）。一些 bin 将包含多个数据点，而其中很多根本不包含任何数据点。
  • 在上面的示例中，bin 的大小/长度（我称之为 windowssize）是 1000。我还不知道最适合我的目的的 windowssize 是多少，所以，我希望能够有一个脚本有一个“可调”的窗口大小。
  • 此外，示例 1-1000, 1001, 2001, [...] 显然没有在 bin 中显示任何重叠。然而，这意味着我失去了很多“敏感性”，因此失去了知识/信息。所以我希望能够创建重叠的窗口/箱。例如，1-1000,501-1500、1501-2500、2001-3000。这些箱有 500 个重叠。我将能够设置箱之间的重叠量。
  • 我想将每个 bin 写入自己的文件，即使 bin 不包含任何数据点。

这是我之前给出的解释，显然不是很好

使用 AWK 我试图在整数列表上“滑动窗口”。如果我拆分这个数据集，每个（可能重叠的）bin 中会有多少个数据点？我喜欢设置 binsize（或 windowssize）并在 bin 之间重叠。这种方法使我能够了解本地数据点的密度。--> 我有一点 AWK 经验，有人告诉我 AWK 应该能够胜任这项工作，我更喜欢使用 AWK。但是，我也对其他想法持开放态度（例如 Python）。

  • 我有一个数据文件，其中包含大约 27,500 个介于 0 和 46041637（4600 万）之间的整数以及每个数据点的描述。
  • 因为我想了解“改变”分辨率的效果，所以我想在不同的窗口大小和各个窗口之间的重叠情况下玩一点。
  • 我喜欢将每个“窗口”的所有内容写入一个单独的文件，并根据“windows-startpoint”命名文件。

我准备了一些示例文件，附在下面。但是，为了更容易理解这里是另一个非常非常简化的示例：

Datarange 1 - 10
Integers in dataset: 2,4,5,6,9
####
####Example 1: Windowsize=5,overlap=2####
####
file name ="1" contents are: (range is 1 - 5)
2
4
5
file name="3" contents are: (range is 4 - 8, that is, two overlap with previous range)
4
5
6
file name="7" contents are: (range is 7 - 10, if the range was larger, it would be 7 - 11)
9
####
####Example 2: Windowsize=3,overlap=0####
####
file name="1" contents are (range 1 - 3)
2
file name="4" contents are (range 4 - 6)
4
5
6
file name="7" contents are (range 7 - 9)
9
file name="9" contents are (range 10 - 10)
<none>

示例输入文件

    3579                    
    3661                    
    3752    EXON    3706    4407    +   Solyc06g005000.2.1.1Solyc06g005000.2.1
    3947    EXON    3706    4407    +   Solyc06g005000.2.1.1Solyc06g005000.2.1
    6734    INTRON  5605    7662    +   Solyc06g005000.2.1.2Solyc06g005000.2.1
    6865    INTRON  5605    7662    +   Solyc06g005000.2.1.2Solyc06g005000.2.1
    6915    INTRON  5605    7662    +   Solyc06g005000.2.1.2Solyc06g005000.2.1
    8961                    
    13471   EXON    13449   13532   +   Solyc06g005020.1.1.2Solyc06g005020.1.1
    13561   INTRON  13533   13710   +   Solyc06g005020.1.1.2Solyc06g005020.1.1
    22226   EXON    22106   22261   +   Solyc06g005030.1.1.1Solyc06g005030.1.1
    22516                   
    22556                   
    36903   INTRON  36836   36915   +   Solyc06g005060.2.1.1Solyc06g005060.2.1
    37377   EXON    36916   37800   +   Solyc06g005060.2.1.2Solyc06g005060.2.1
    37605   EXON    36916   37800   +   Solyc06g005060.2.1.2Solyc06g005060.2.1
    37935   3P_UTR  37801   38132   +   Solyc06g005060.2.1.0Solyc06g005060.2.1
    167942  5P_UTR  167930  167956  -   Solyc06g005140.2.1.0Solyc06g005140.2.1
    168020  INTRON  167957  169025  -   Solyc06g005140.2.1.2Solyc06g005140.2.1
    168153  INTRON  167957  169025  -   Solyc06g005140.2.1.2Solyc06g005140.2.1

具有不同窗口大小和重叠的示例输出文件

> AWK -v windowsize=50000 -v overlap=0 -f awkscript input.file
> ls
1      50001
100001 150001
> cat 1
    3579                    
    3661                                        
    3752    EXON    3706    4407    +   Solyc06g005000.2.1.1Solyc06g005000.2.1
    3947    EXON    3706    4407    +   Solyc06g005000.2.1.1Solyc06g005000.2.1
    6734    INTRON  5605    7662    +   Solyc06g005000.2.1.2Solyc06g005000.2.1
    6865    INTRON  5605    7662    +   Solyc06g005000.2.1.2Solyc06g005000.2.1
    6915    INTRON  5605    7662    +   Solyc06g005000.2.1.2Solyc06g005000.2.1
    8961                    
    13471   EXON    13449   13532   +   Solyc06g005020.1.1.2Solyc06g005020.1.1
    13561   INTRON  13533   13710   +   Solyc06g005020.1.1.2Solyc06g005020.1.1
    22226   EXON    22106   22261   +   Solyc06g005030.1.1.1Solyc06g005030.1.1
    22516                   
    22556                   
    36903   INTRON  36836   36915   +   Solyc06g005060.2.1.1Solyc06g005060.2.1
    37377   EXON    36916   37800   +   Solyc06g005060.2.1.2Solyc06g005060.2.1
    37605   EXON    36916   37800   +   Solyc06g005060.2.1.2Solyc06g005060.2.1
    37935   3P_UTR  37801   38132   +   Solyc06g005060.2.1.0Solyc06g005060.2.1                  
    3752    EXON    3706    4407    +   Solyc06g005000.2.1.1Solyc06g005000.2.1
    3947    EXON    3706    4407    +   Solyc06g005000.2.1.1Solyc06g005000.2.1
    6734    INTRON  5605    7662    +   Solyc06g005000.2.1.2Solyc06g005000.2.1
    6865    INTRON  5605    7662    +   Solyc06g005000.2.1.2Solyc06g005000.2.1
    6915    INTRON  5605    7662    +   Solyc06g005000.2.1.2Solyc06g005000.2.1
    8961                    
    13471   EXON    13449   13532   +   Solyc06g005020.1.1.2Solyc06g005020.1.1
    13561   INTRON  13533   13710   +   Solyc06g005020.1.1.2Solyc06g005020.1.1
    22226   EXON    22106   22261   +   Solyc06g005030.1.1.1Solyc06g005030.1.1
    22516                   
    22556                   
    36903   INTRON  36836   36915   +   Solyc06g005060.2.1.1Solyc06g005060.2.1
    37377   EXON    36916   37800   +   Solyc06g005060.2.1.2Solyc06g005060.2.1
    37605   EXON    36916   37800   +   Solyc06g005060.2.1.2Solyc06g005060.2.1
    37935   3P_UTR  37801   38132   +   Solyc06g005060.2.1.0Solyc06g005060.2.1
 > cat 50001

 > cat 100001

 > cat 150001
    167942  5P_UTR  167930  167956  -   Solyc06g005140.2.1.0Solyc06g005140.2.1
    168020  INTRON  167957  169025  -   Solyc06g005140.2.1.2Solyc06g005140.2.1
    168153  INTRON  167957  169025  -   Solyc06g005140.2.1.2Solyc06g005140.2.1

> #And with some different paramenters
> AWK -v windowsize=160000 -v overlap=10000 -f awkscript input.file
> ls
1    10001
> cat 1
    3579                    
    3661                    
    3752    EXON    3706    4407    +   Solyc06g005000.2.1.1Solyc06g005000.2.1
    3947    EXON    3706    4407    +   Solyc06g005000.2.1.1Solyc06g005000.2.1
    6734    INTRON  5605    7662    +   Solyc06g005000.2.1.2Solyc06g005000.2.1
    6865    INTRON  5605    7662    +   Solyc06g005000.2.1.2Solyc06g005000.2.1
    6915    INTRON  5605    7662    +   Solyc06g005000.2.1.2Solyc06g005000.2.1
    8961                    
    13471   EXON    13449   13532   +   Solyc06g005020.1.1.2Solyc06g005020.1.1
    13561   INTRON  13533   13710   +   Solyc06g005020.1.1.2Solyc06g005020.1.1
    22226   EXON    22106   22261   +   Solyc06g005030.1.1.1Solyc06g005030.1.1
    22516                   
    22556                   
    36903   INTRON  36836   36915   +   Solyc06g005060.2.1.1Solyc06g005060.2.1
    37377   EXON    36916   37800   +   Solyc06g005060.2.1.2Solyc06g005060.2.1
    37605   EXON    36916   37800   +   Solyc06g005060.2.1.2Solyc06g005060.2.1
    37935   3P_UTR  37801   38132   +   Solyc06g005060.2.1.0Solyc06g005060.2.1
> cat 10001
    13471   EXON    13449   13532   +   Solyc06g005020.1.1.2Solyc06g005020.1.1
    13561   INTRON  13533   13710   +   Solyc06g005020.1.1.2Solyc06g005020.1.1
    22226   EXON    22106   22261   +   Solyc06g005030.1.1.1Solyc06g005030.1.1
    22516                   
    22556                   
    36903   INTRON  36836   36915   +   Solyc06g005060.2.1.1Solyc06g005060.2.1
    37377   EXON    36916   37800   +   Solyc06g005060.2.1.2Solyc06g005060.2.1
    37605   EXON    36916   37800   +   Solyc06g005060.2.1.2Solyc06g005060.2.1
    37935   3P_UTR  37801   38132   +   Solyc06g005060.2.1.0Solyc06g005060.2.1
    167942  5P_UTR  167930  167956  -   Solyc06g005140.2.1.0Solyc06g005140.2.1
    168020  INTRON  167957  169025  -   Solyc06g005140.2.1.2Solyc06g005140.2.1
    168153  INTRON  167957  169025  -   Solyc06g005140.2.1.2Solyc06g005140.2.1

非常感谢您的帮助！

我最初的问题几乎没有调整，因为我需要的计算时间比我预期的要多。

是否有可能不是将特定窗口中的所有记录写入其自己的文件，而是将每个窗口的“统计信息”写入表中的一行？对于统计数据，我的意思是特定窗口包含多少条记录以及每种类型有多少条记录。应用到上面的示例中，它看起来像这样：

> python script.py 160000 10000 file (using the script from sidharth c nadhan)
> cat result
window | total | exons | intron | 3P_UTR | 5P_UTR
1      |  17   |   6   |   5    |    1   |   0
10001  |  12   |   4   |   4    |    1   |   1
4

1 回答 1

1

尝试这个 ：

import sys,os,collections
list1,set1=list(),set()
dict1 = collections.defaultdict(list)
dict2 = collections.defaultdict(int)
wind , overl, maxim= int(sys.argv[1]),int(sys.argv[2]),int(sys.argv[4])
for line in open(sys.argv[3]):
  try : set1.add(line.split()[1])
  except : pass
  for i in xrange(1,maxim,wind-overl):
    if int(line.split()[0]) in xrange(i,wind+i): dict1[i].append(line)
print "Window\tTotal\t",
for entri in set1 : print entri,"\t",
print "\n",
for j in xrange(1,maxim,wind-overl):
  dict2.clear()
  for line in dict1[j]:
    try : dict2[line.split()[1]] +=1
    except: pass
  print j,"\t",len(dict1[j]),"\t",
  for entri in set1:
    print dict2[entri],"\t",
  print "\n",

用法：python script.py窗口重叠文件 maxim

其中 maxim 是指输入文件中的最大整数。示例输入文件中的 maxim = 168153。将其作为命令行参数提供可以提高计算速度。

于 2013-05-20T11:07:21.480 回答