从TCGA数据中提取lncRNA的表达量时,需要知道lncRNA的编号和对应的名称。这些信息可以从GTF文件中提取。
提取的话,可以采用如下的代码实现。
#!/usr/bin/perl -w use strict; my $biotype_file = shift @ARGV; my $gtf = shift @ARGV; my $biotype = shift @ARGV; my %biotype_list; open my $fh1, $biotype_file or die; while (<$fh1>) { chomp; my @array = split /\t/, $_; if($array[2]eq $biotype){ $biotype_list{$array[0]} = 1; } } close $fh1; open my $out, ">${biotype}_info.txt" or die; print $out "Gene_id\tGene_id_info\tgene_name\tbiotype\n"; open my $fh2, $gtf or die; while (<$fh2>) { chomp; next if /^#/; my @array = split /\t/, $_; next unless ($array[2] eq "gene"); $array[8] =~ /gene_id\s+"(\S+?)";.*gene_type\s+"(\S+?)";.*gene_name\s+"(\S+?)";/; my $geneid = $1; my $genebiotype = $2; my $genename = $3; my $gene_id_norm=(split("\\.",$geneid))[0]; if ($biotype_list{$genebiotype}) { print $out "$gene_id_norm\t$geneid\t$genename\t$genebiotype\n"; } } close $fh2;
如果您对TCGA数据挖掘感兴趣,请学习我们的TCGA相关课程:
《GSEA富集分析》
如果觉得我的文章对您有用,请随意打赏。你的支持将鼓励我继续创作!