- 論壇徽章:
- 0
|
本帖最后由 弦斷有誰聽1053476508 于 2014-12-15 16:53 編輯
任務(wù):步驟1,Orthologs-msu.txt是原始文件,讓其每一行行生成1個ID文件(如ORTHOMCL2.txt) ,文件名是每一行第一個元素。
步驟2,然后以文件名為散列的key,文件內(nèi)容為元素構(gòu)造散列保存。
步驟3,再用每個ID文件的每一行(如下面的APK_ORTHOMCL2),從同一個Allseq.fa(DNA文件)中找到對應(yīng)的APK_ORTHOMCL2提取DNA序列
如,APK_ORTHOMCL2
>AT1G51370.2 | Symbols: | F-box/RNI-like/FBD-like domains-containing protein | chr1:19045615-19046825 FORWARD LENGTH=1118
ATGGTGGGTGGCAAGAAGAAAACCAAGATATGTGACAAAGTGTCACATGAGGAAGATAGGATAAGCCAGTTACCGGAACC
TTTGATATCTGAAATACTTTTTCATCTTTCTACCAAGGACTCTGTCAGAACAAGCGCTTTGTCTACCAAATGGAGATATC
TTTGGCAATCGGTTCCTGGATTGGACTTAGACCCCTACGCATCCTCAAATACCAATACAATTGTGAGTTTTGTTGAAAGT
TTTTTTGATTCCCACAGGGATTCATGGATACGCAAACTCCGTTTAGATTTGGGTTATCATCATGATAAGTATGATCTCAT
GTCATGGATTGATGCTGCGACTACGCGTAGGATTCAGCATCTTGATGTTCATTGTTTTCACGATAATAAGATACCCTTGA
GCATATATACATGCACGACGTTGGTACACTTACGACTCCGTTGGGCTGTCTTGACTAATCCCGAGTTTGTTTCCTTACCT
生成N個fa序列文件,找不到對應(yīng)ID的生產(chǎn)空文件。
原始文件就像下面的內(nèi)容:
orthologous_group_ID num_orthologs/in-paralogs num_species species orthologs/in-paralogs
APK_ORTHOMCL0 416 2 rice sorghum LOC_Os04g04140 LOC_Os05g29160 LOC_Os10g02640
APK_ORTHOMCL2 212 3 brachy rice sorghum Bradi1g21530 Bradi2g26540 Bradi2g49480 Bradi2g61260
其中要求生成的ID文件如下:
APK_ORTHOMCL0
416
2
rice sorghum
LOC_Os04g04140
LOC_Os05g29160
LOC_Os10g02640
別人寫的步驟一的一段代碼,錯誤出在哪,怎么改?或者高手給重寫一段。
open (IN,"<$ARGV[0]");
my ($i,@in,$linename,@id);
while (<IN>)
{
@in = split("\s", $_);
$linename = shift (@in);
@id = join ("\n", @in);
print "$linename\n@id";
$i++;
open my $out, ">$i";
print $out sth
步驟二沒有參考腳本
步驟三的參考腳本如下:
#tempo 從fasta格式中用名稱提取序列
#cd C:\strawberry\perl\bin
# perl ID-SEQok.pl ORTHOMCL2.txt zmossbat.fa >mcl2.fa
#!/usr/bin/perl -w
#use strict;
#use warnings;
($ARGV[0] && $ARGV[1]) or die "usage: $0 NameFile DataBase\nShortName to FullName\n";
my $nameFile=$ARGV[0];
my $dataFile=$ARGV[1];
open(my $nf,$nameFile) or die "cannot open $nameFile";
open(my $df,$dataFile) or die "cannot open $dataFile";
undef $/;
my $temp;
$temp=<$nf>;
my @names=sort split /\n/,$temp;
close($nf);
foreach my $name (@names)
{
$name=~s/\|/\\|/g;
}
$/=">";
while($temp=<$df>)
{
chomp($temp);
$temp=~m/(.*)/;
my $head=$1;
foreach my $name(@names)
{
if($head=~m/$name/i)
{
print ">$temp";
}
}
}
close($df);
我是植物學(xué)碩士,但是對計算機語言沒什么造詣,高手幫幫忙吧,最好把三個步驟連在一起寫一個程序。
---------------------------------------------------------------------------------------------------------------------------------------------------------
12月15編輯:
原始文件(Orthologs-msu.txt有約2000行)和要求生成約2000個fa序列文件,過程如
140405lhnzm2uvtzp522vm.jpg (803.31 KB, 下載次數(shù): 65)
下載附件
2014-12-15 16:52 上傳
Allseq.fa文件有200g,太大沒法上傳,只能發(fā)測試文件。測試文件在下面的附件中,具體測試就是使原始文件Orthologs-msu.txt每一行生成1個ID文件,文件名是每一行第一個元素,文件內(nèi)容為ID。ID文件的每一行為一個子ID(如ORTHOMCL2.txt),再用每個ID文件的每一行(如下面的LOC_Os04g04140),從同一個Allseq.fa(DNA文件)中找到對應(yīng)的(如LOC_Os04g04140)提取序列,生成N個fa序列文件。效果如APK_ORTHOMCL0.fa和APK_ORTHOMCL2.fa,提取不到對應(yīng)序列的如APK_ORTHOMCL1.fa為空!
測試文件(12月15).zip
(17.53 KB, 下載次數(shù): 4)
2014-12-15 14:01 上傳
點擊文件名下載附件
|
|