add intron/exon schema and mrna cds loader

[dbsnp.git] / utils / human_mrna_cds_insert.pl

diff --git a/utils/human_mrna_cds_insert.pl b/utils/human_mrna_cds_insert.pl
new file mode 100755 (executable)
index 0000000..ad15784
--- /dev/null
+++ b/utils/human_mrna_cds_insert.pl
@@ -0,0 +1,150 @@
+#!/usr/bin/perl
+
+use warnings;
+use strict;
+
+# psql snp -c "COPY mrna_cds_table FROM STDIN WITH DELIMITER '   ' NULL AS 'NULL'";
+
+my $prot_table;
+
+while (<DATA>) {
+    chomp;
+    my ($codon,$aa,$aa_med,$aa_long) = split /\t/;
+    $prot_table->{$codon} = $aa;
+}
+
+my $current_stanza;
+while (<>) {
+    $current_stanza .= $_;
+    if ($_ =~ m{^//$}) {
+       handle_mrna_stanza($current_stanza);
+       $current_stanza = '';
+    }
+}
+
+sub handle_mrna_stanza{
+    my ($current_stanza) = @_;
+
+    my ($accession,$version,$gi) =
+       $current_stanza =~ m/^VERSION\s+([^\.]+)\.(\d+)\s+GI:(\d+)$/m;
+    my ($gene) =
+       $current_stanza =~ m/^\s+\/gene=\"?([^"]+)\"?$/m;
+    my ($cds) =
+       $current_stanza =~ m/^\s+CDS\s+(\d+(?:\.\.\d+))$/m;
+    return if not defined $cds;
+    my ($start,$stop) = split /\.\./,$cds;
+    return if not defined $start;
+    if (not defined $stop) {
+       $stop = $start;
+    }
+    my ($sequence) =
+       $current_stanza =~ m/^(ORIGIN.+)$/sm;
+    $sequence =~ s{\s*//\s*$}{};
+    $sequence =~ s{^ORIGIN}{};
+    $sequence =~ s{[\d\s\n]+}{}g;
+    $sequence = uc($sequence);
+    my $sequence2 = substr $sequence,($start-1),($stop-$start+1);
+    # There are other posibilities, like CTG, ACG, etc; so we can't test this.
+    # if ($sequence2 !~ /^A[UT]G/i) {
+    #  print STDERR $sequence2."\n";
+    #  print STDERR $gene."\n";
+    #  die "Sequence doesn't start with a start codon";
+    # }
+    if ($sequence2 !~ /(?:[TU]AA|[TU]AG|[TU]GA)$/) {
+       print STDERR $sequence2."\n";
+       print STDERR $gene."\n";
+       die "Sequence doesn't end with a stop codon";
+    }
+    # figure out cds gi
+    my $cds_gi='NULL';
+    my ($cds_stanza) = $current_stanza =~ /^\s{5}CDS\s+(.+?)^(?:\s{5}\w+|ORIGIN)/ms;
+    if (length $cds_stanza) {
+       ($cds_gi) = $cds_stanza =~ /db_xref="GI:(\d+)"/;
+       if (not defined $cds_gi) {
+           $cds_gi = 'NULL';
+       }
+    }
+    print join("\t",$gi,$cds_gi,$accession,$version,$gene,$start,$stop,length($sequence2), uc($sequence2),convert_to_prot_seq($sequence2))."\n";
+}
+
+sub convert_to_prot_seq{
+    my ($mrna) = @_;
+$mrna =~ s/U/T/g;
+    my $prot = '';
+    for my $pos (0..(length($mrna) / 3 - 1)) {
+       my $codon = substr($mrna,$pos*3,3);
+       if (not exists $prot_table->{$codon}) {
+          $prot .= 'N';
+       }
+       else {
+           $prot .= $prot_table->{$codon};
+       }
+    }
+    return $prot;
+}
+
+__DATA__
+ATT    I       Ile     Isoleucine
+ATC    I       Ile     Isoleucine
+ATA    I       Ile     Isoleucine
+CTT    L       Leu     Leucine
+CTC    L       Leu     Leucine
+CTA    L       Leu     Leucine
+CTG    L       Leu     Leucine
+TTA    L       Leu     Leucine
+TTG    L       Leu     Leucine
+GTT    V       Val     Valine
+GTC    V       Val     Valine
+GTA    V       Val     Valine
+GTG    V       Val     Valine
+TTT    F       Phe     Phenylalanine
+TTC    F       Phe     Phenylalanine
+ATG    M       Met     Methionine
+TGT    C       Cys     Cysteine
+TGC    C       Cys     Cysteine
+GCT    A       Ala     Alanine
+GCC    A       Ala     Alanine
+GCA    A       Ala     Alanine
+GCG    A       Ala     Alanine
+GGT    G       Gly     Glycine
+GGC    G       Gly     Glycine
+GGA    G       Gly     Glycine
+GGG    G       Gly     Glycine
+CCT    P       Pro     Proline
+CCC    P       Pro     Proline
+CCA    P       Pro     Proline
+CCG    P       Pro     Proline
+ACT    T       Thr     Threonine
+ACC    T       Thr     Threonine
+ACA    T       Thr     Threonine
+ACG    T       Thr     Threonine
+TCT    S       Ser     Serine
+TCC    S       Ser     Serine
+TCA    S       Ser     Serine
+TCG    S       Ser     Serine
+AGT    S       Ser     Serine
+AGC    S       Ser     Serine
+TAT    Y       Tyr     Tyrosine
+TAC    Y       Tyr     Tyrosine
+TGG    W       Trp     Tryptophan
+CAA    Q       Gln     Glutamine
+CAG    Q       Gln     Glutamine
+AAT    N       Asn     Asparagine
+AAC    N       Asn     Asparagine
+CAT    H       His     Histidine
+CAC    H       His     Histidine
+GAA    E       Glu     Glutamic acid
+GAG    E       Glu     Glutamic acid
+GAT    D       Asp     Aspartic acid
+GAC    D       Asp     Aspartic acid
+AAA    K       Lys     Lysine
+AAG    K       Lys     Lysine
+CGT    R       Arg     Arginine
+CGC    R       Arg     Arginine
+CGA    R       Arg     Arginine
+CGG    R       Arg     Arginine
+AGA    R       Arg     Arginine
+AGG    R       Arg     Arginine
+TAA    *       *       Stop codon
+TAG    *       *       Stop codon
+TGA    *       *       Stop codon