f*****a
发帖数: 156 | 1 编程水平不行，用Python的readline和C++的getline写的code慢死了，求高效算法。
问题如下：
两个输入文件，文件1是fastq格式，几千万行，每行最多100个字符
（每4行是一个序列的信息，其中的第2行是序列本身）:
@SEQ_ID1
GATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
+
!''*((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCC65
@SEQ_ID2
AATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
+
+''*((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCC65
@SEQ_ID3
TAGGCGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
+
C''F((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCCAA
...
文件2是两列的文本，1000行：
BC1 GATTT
BC2 TAGGC
....
程序1：对于文件2中每个BC，计算文件1中有多少序列是以该BC的序列开头的。
例如：对于BC1，在文件1中的第4N+2行中有多少是GATTT开头的。
程序2：和程序1有点类似，不过要求是对于文件2中每个BC，
把文件1中以该BC序列开头的所有序列的4行信息写入单独的文件。
例如：对于BC1，在文件1中的第4N+2行中查找所有是GATTT开头的
把每个匹配的4行信息写入文件BC1.txt
第一个匹配写入的是：
@SEQ_ID1
GATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
+
!''*((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCC65 |
w****w
发帖数: 521 | |
f*****a
发帖数: 156 | 3 不是的，从6到15不等。
我现在的方法是对文件1循环，每次读4行。同时建立1000个BC#.txt文件，比较每个BC
和读入的第二行，再计数或写入。程序可以运行，不过很慢。
[发表自未名空间手机版 - m.mitbbs.com]
【在 w****w 的大作中提到】
: 文件2中第二列长度都是6？
|
w****w
发帖数: 521 | 4 最好把它们都搞成15，这样只要查一次就可以了。
以第二列为key, (count,file)为value建dictionary,然后取seq开始的6，..., 15，查
dictionary。
BC
【在 f*****a 的大作中提到】
: 不是的，从6到15不等。
: 我现在的方法是对文件1循环，每次读4行。同时建立1000个BC#.txt文件，比较每个BC
: 和读入的第二行，再计数或写入。程序可以运行，不过很慢。
:
: [发表自未名空间手机版 - m.mitbbs.com]
|
f*****a
发帖数: 156 | 5 不太明白，可否给个python的code? 多谢!
[发表自未名空间手机版 - m.mitbbs.com]
【在 w****w 的大作中提到】
: 最好把它们都搞成15，这样只要查一次就可以了。
: 以第二列为key, (count,file)为value建dictionary,然后取seq开始的6，..., 15，查
: dictionary。
:
: BC
|
k**********g
发帖数: 989 | 6 http://stackoverflow.com/questions/843154/fastest-way-to-find-t
As a sidenote,
(1) Make sure your text file is not Unicode, and make sure its end-of-line
marker is consistent. The processing of variable-length Unicode-encoded
character is what prevents the standard library writers from making
optimizations.
(2) If the BC sequences have very short length (such as 5), and the
alphabets allowed in BC sequences are just {C, G, T, A} then the exhaustive
number of possible BC sequences is only 4 to the power of 5 = 1024. Thus
you can build a histogram of all BC sequences (whether they appear in File 2
or not) with length 5 or less. |
w****w
发帖数: 521 | 7 http://webbew.dyx.com/forum/viewtopic.php?f=2&t=8
【在 f*****a 的大作中提到】
: 不太明白，可否给个python的code? 多谢!
:
: [发表自未名空间手机版 - m.mitbbs.com]
|
m*****n
发帖数: 204 | 8 If you build a prefix tree out of the BCs first, you'll only need to walk
the tree once to find the match.
BC
【在 f*****a 的大作中提到】
: 不是的，从6到15不等。
: 我现在的方法是对文件1循环，每次读4行。同时建立1000个BC#.txt文件，比较每个BC
: 和读入的第二行，再计数或写入。程序可以运行，不过很慢。
:
: [发表自未名空间手机版 - m.mitbbs.com]
|
|
d****n
发帖数: 1637 | 9 FASTQ file customize deindex?
step1, use klib->kseq to read fastq files. No other lib is faster than this.
step 2, build a hash base on your barcode list, and open 1000 file handlers
as the hash values.
step3. iterate through your fastq files,(using kseq), when a
barcode match a key write it to the value(file hanlder)
I bet you no other solution is faster than this.
【在 f*****a 的大作中提到】
: 编程水平不行，用Python的readline和C++的getline写的code慢死了，求高效算法。
: 问题如下：
: 两个输入文件，文件1是fastq格式，几千万行，每行最多100个字符
: （每4行是一个序列的信息，其中的第2行是序列本身）:
: @SEQ_ID1
: GATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
: +
: !''*((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCC65
: @SEQ_ID2
: AATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
|
o*******p
发帖数: 27 | 10 如果长度不等，可以用10 (=15-6+1)个dictionary，每个对应相对的长度：
ds = []
for i in range(10): ds.append({})
for line in open('file2.txt').readlines():
name, seq = line.strip().split()
length = len(seq)
ds[length - 6][seq] = 0
这样第一个文件的每一个序列需要与ds的每一个dictionary相比，对应不同的长度：
for seq in seqs: #省略了FASTQ reading
for i in range(10):
prefix = seq[6 + i]
if prefix in ds[i]:
ds[i] += 0
这样比直接比较应该快很多。
【在 w****w 的大作中提到】
: 最好把它们都搞成15，这样只要查一次就可以了。
: 以第二列为key, (count,file)为value建dictionary,然后取seq开始的6，..., 15，查
: dictionary。
:
: BC
|
|
|
f*****a
发帖数: 156 | |
f*****a
发帖数: 156 | 12 多谢。没用过klib。
这个和前面的python算法基本一样
如果用c++来算会更快么？
this.
handlers
【在 d****n 的大作中提到】
: FASTQ file customize deindex?
: step1, use klib->kseq to read fastq files. No other lib is faster than this.
: step 2, build a hash base on your barcode list, and open 1000 file handlers
: as the hash values.
: step3. iterate through your fastq files,(using kseq), when a
: barcode match a key write it to the value(file hanlder)
: I bet you no other solution is faster than this.
|
a****e
发帖数: 9589 | 13 用shell script 就可以搞定，还快
【在 f*****a 的大作中提到】
: 编程水平不行，用Python的readline和C++的getline写的code慢死了，求高效算法。
: 问题如下：
: 两个输入文件，文件1是fastq格式，几千万行，每行最多100个字符
: （每4行是一个序列的信息，其中的第2行是序列本身）:
: @SEQ_ID1
: GATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
: +
: !''*((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCC65
: @SEQ_ID2
: AATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
|
b***i
发帖数: 3043 | 14 是否可以读入到内存中？
是否可以使用Trie?
【在 f*****a 的大作中提到】
: 编程水平不行，用Python的readline和C++的getline写的code慢死了，求高效算法。
: 问题如下：
: 两个输入文件，文件1是fastq格式，几千万行，每行最多100个字符
: （每4行是一个序列的信息，其中的第2行是序列本身）:
: @SEQ_ID1
: GATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
: +
: !''*((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCC65
: @SEQ_ID2
: AATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
|
d******e
发帖数: 2265 | 15 这个
shang sed 很简单吧
【在 f*****a 的大作中提到】
: 编程水平不行，用Python的readline和C++的getline写的code慢死了，求高效算法。
: 问题如下：
: 两个输入文件，文件1是fastq格式，几千万行，每行最多100个字符
: （每4行是一个序列的信息，其中的第2行是序列本身）:
: @SEQ_ID1
: GATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
: +
: !''*((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCC65
: @SEQ_ID2
: AATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
|
b*****e
发帖数: 474 | 16 才几千万行，直接排序再分段分成很多小文件不就行了吗。
只需干一次，以后就很快了 |
