Re^2: Iso electric point calculation using perl

[yuvraj_ghaly]

See As the caption of my talk is to calculate Iso electric point of protein sequences which you havenot understood properly.

I have a Fasta file containing numerous protein sequences in Fasta format. You can check NCBI homepage for that. www.ncbi.nlm.nih.gov. My task is to calculate Isoelectric point of each protein sequences. The script I composed is not helping me. I'm not a Perl expert so I commited some mistake. I really dont know what. I keep on trying but I'm keep messing around. I need expert advice.

[Anonymous Monk]

repost what you have after corrections.What you posted the first time makes no sense.

[yuvraj_ghaly]

 may be this code help you to know what I'm trying to do</p>
<code>
#! /usr/bin/perl -w
use strict;
$ph=0;
$qn1='';
$qn2='';
$qn3='';$qn4='';$qn5='';
$qp1='';$qp2='';$qp3='';$qp4='';
$nq='';
open (S, "$ARGV[0]") || die "cannot open FASTA file to read: $!";

my %s;# a hash of arrays, to hold each line of sequence
my %seq; #a hash to hold the AA sequences.
my $key;

while (<S>){ #Read the FASTA file.
    chomp;

    if (/>/){
        s/>//;
        $key= $_;
    }else{
        push (@{$s{$key}}, $_);
    }

}

foreach my $a (keys %s){
    my $s= join("", @{$s{$a}});
    $seq{$a}=$s;
    #print("$a\t$s\n");
}

my @aa= qw(Length);
#my @aa= qw(A R N D C Q E G H I L K M F P S T W Y V);
my $aa= join("\t", @aa);
print ("Sequence\t$aa\n");

foreach my $k (keys %seq){
    my %count; # a hash to hold the count for each amino acid in the p
+rotein.
    my @seq= split(//, $seq{$k});
    foreach my $r(@seq){
        $count{$r}++;
    }
my @row;
push(@row, $k);
    foreach my $a (@aa){
        $count{$a}||=0;
        $count{$a}=length($seq{$k});
 #       $count{$a}= sprintf("%0.6f",($count{$a}/length($seq{$k})));
        push(@row,$count{$a});
    }
my $row= join("\t",@row);
print("\n$row\n");
}

my $pkaC= 9.0; 
my $pkaD= 4.05; 
my $pkaE= 4.45; 
my $pkaY= 10.0;
my $pkaH= 5.98;
my $pkaK= 10.0;
my $pkaR= 12.0;
foreach my $pka (@row){
$qn1= -1/(1+10^(3.65-$ph));
$qn2= -$cys/(1+10^(9.0-$ph));
$qn3= -$asp/(1+10^(4.05-$ph));
$qn4= -$glu/(1+10^(4.45-$ph));
$qn5= -$tyr/(1+10^(10.0-$ph));
$qp1= 1/(1+10^($ph-8.2));
$qp2= $his/(1+10^($ph-5.98));
$qp3= $lys/(1+10^($ph-10.0));
$qp4= $arg/(1+10^($ph-12.0));
$nq= $qn1+$qn2+$qn3+$qn4+$qn5+$qp1+$qp2+$qp3+$qp4;
}

if ($ph>=7){
next;
}
if ($nq<=0){
next;
$ph++= 0.01;
}
open (OUT, ">"."$nq"."_"."$number".".xls") or die "Cannot open file!";
print OUT "$nq";

[Anonymous Monk]

What does this mean?

@aa= qw(Length);

Because it's not ok.

[Anonymous Monk]

correct this:

my @aa= qw(A R N D C Q E G H I L K M F P S T W Y V);
my $aa= join("\t", @aa);

assign a new variable for $aa.

and this:

$ph++= 0.01;

into:

$ph+= 0.01;

[Anonymous Monk]

and plz. explain what it is you're trying to do here:

$qn2= -$cys/(1+10^(9.0-$ph));

[Anonymous Monk]

And this needs explanation:

while (<S>){ #Read the FASTA file.
    chomp;

    if (/>/){
        s/>//;
        $key= $_;
    }else{
        push (@{$s{$key}}, $_);
    }

}

[Anonymous Monk]

And this needs explanation:

while (<S>){ #Read the FASTA file.
    chomp;

    if (/>/){
        s/>//;
        $key= $_;
    }else{
        push (@{$s{$key}}, $_);
    }

}

[Anonymous Monk]

After having a good look at your entire program i must conclude i cannot help you unless i learn biology myself. good luck.