Sequencing Analysis - Problem Verarbeitung der Masse (Allgemeines zu Perl)

[thread]17018[/thread]

Sequencing Analysis - Problem Verarbeitung der Masse

Tags: perl5 Ähnliche Threads

Leser: 18

Articles: hide open all | hide show old branches

+4 replies

Bells

2011-12-26 17:24

User since
2011-12-26
2 Artikel
BenutzerIn
[default_avatar]

So, ich bin ehrlich, dass ich so meine Probleme mit dem Programmieren habe. Wir haben ein Programmierprojekt an der Uni, für Sequencing Analysis. Die Aufgabe bestand darin: Wir haben Reads und unterschiedliche Chromosomen bekommen. Alles fein und schick, die soll man gegeneiander alignen, dafür haben wir in der Gruppe den Needlerman-Wunsch Algorithmus zu einem Programm umgeschrieben, was auch so läuft das es alignen tut und könnte in der Theorie nach beim raussuchen behilflich werden.Es kann den Match Score berechen, denGap Scoreund den Best Score raussuchen kann.Alles schön und gut - jetzt setzt das Problem ein.

Erst einmal ist die Masse an Daten viel zu groß - das Programm kann bei Eingabe im Terminal nur begrenzt Daten aufnehmen. Wir sollen mehere Chromosomen untersuchen, die reads die zum absuchen gedacht sind, sind 4000 Basen allein bereits lang.

Die Masse der Daten,das Programm kann nur eine begrenzte Daten überhaupt berechnen und gegenseitig alginen.Das ist schon einmal ein Problem. Mir stellt sich die Frage, ob man ein Programm überhaupt so umstricken kann das es auch die Masse an Daten abarbeiten tut und die Ergebnisse dann in einer weiteren Datei abspeichern tut. Denn zur Aufgabe gehört es, wir sollen eine Datei erstellen, Format ist uns überlassen, wo die Ergebnisse anschaulich dargestellt sind.

Alle Datein, mit den Basen und Sequenzen sind im Fastaformat vorgegeben. Ich bin ehrlich, wir sind wirkliche Perl - Noobs, und habe zu zweit an dem Programm die 3 Tage über Weihnachten gehockt und sind froh, dass es so erstmal läuft, aber das klappt einfach vorn hinten nicht mit dem abarbeiten der Masse an Daten die es zu bearbeiten gibt und irgedwie müsste das Programm rund um die Uhr laufen und die Abgleiche machen und irgendwie speichern in einer Datei. Nur - wie soll das gehen? Da fehlt jeder Ansatz.
Last edited: 2011-12-26 19:05:12 +0100 (CET)

+3 replies

Bells

2011-12-26 20:00

User since
2011-12-26
2 Artikel
BenutzerIn
[default_avatar]

Das wäre das Programm - im Moment halt mit Eingabe übers Terminal. Vielleicht wäre es sinnig - Fasta.datei mit den beiden Sequenz auslesen lassen.. damit das Programm dann arbeiten kann und dann in einer neuen Datei dann die Ergenisse liefert?! Hier ist die Frage wie setzt man das Codetechnisch um? :(

und von den Modulen habe ich gelesen Bio::Fasta etc. - möchte der jeweilige Dozent nicht haben. Und wir haben jetzt zwar einiges gelesen dazu, aber ein Ansatz fehlt mir jedoch wie ich wo anfange.

Das dies:

Code (perl): (dl )

print"Gib erste Sequenenz ein \n";
my $seq1=<STDIN>;

print"Gib  zweite Sequenz ein \n";
my $seq2=<STDIN>;

Darüber bin ich mir im klaren,weil so nur kleine Sequezen per Hand eingegeben werden..und der Wusch ist ja die Reads über die chromosomen laufezu lassen,um durch Matching..das passende Chromosom herausfinden zu können.

Code (perl): (dl )

use strict;
use warnings;
#debug = 1 displays score and traceback matrix


my $debug = 0;
my $now = time;

print"Gib erste Sequenenz ein \n";
my $seq1=<STDIN>;

print"Gib  zweite Sequenz ein \n";
my $seq2=<STDIN>;

#################################################################################################### 
#                                                       Schema für Score                                                                                                          #
####################################################################################################

my $MATCH    =  1; # +1 für Buchtstabe, der ein Match macht
my $MISMATCH = -1; # -1 für Buchstabe, der einen Mismatch macht
my $GAP      = -1; # -1 für Lücke
 
 
#################################################################################################### 
#                      Initialisierung / Programmstart                                             #
####################################################################################################
my @matrix;
$matrix[0][0]{score}   = 0;
$matrix[0][0]{pointer} = "none";
for(my $j = 1; $j <= length($seq1); $j++) {
$matrix[0][$j]{score}   = $GAP * $j;
$matrix[0][$j]{pointer} = "left";
}
for (my $i = 1; $i <= length($seq2); $i++) {
$matrix[$i][0]{score}   = $GAP * $i;
$matrix[$i][0]{pointer} = "up";
}

####################################################################################################
#                                                                       Auffüllen                                                      #
####################################################################################################
for(my $i = 1; $i <= length($seq2); $i++) {
for(my $j = 1; $j <= length($seq1); $j++) {
my ($diagonal_score, $left_score, $up_score);
 
 
####################################################################################################
#                            berechnet match score                                                                                                 #
####################################################################################################
my $letter1 = substr($seq1, $j-1, 1);
my $letter2 = substr($seq2, $i-1, 1);
if ($letter1 eq $letter2) {
$diagonal_score = $matrix[$i-1][$j-1]{score} + $MATCH;
}
else {
$diagonal_score = $matrix[$i-1][$j-1]{score} + $MISMATCH;
}
 
 
####################################################################################################
#                                                       berechnet gap scores                                                                                               #
####################################################################################################
$up_score   = $matrix[$i-1][$j]{score} + $GAP;
$left_score = $matrix[$i][$j-1]{score} + $GAP;



#################################################################################################### 
#                             aussuchen des best score                                                                                     #
####################################################################################################
if ($diagonal_score >= $up_score) {
if ($diagonal_score >= $left_score) {
$matrix[$i][$j]{score}   = $diagonal_score;
$matrix[$i][$j]{pointer} = "diagonal";
}
else {
$matrix[$i][$j]{score}   = $left_score;
$matrix[$i][$j]{pointer} = "left";
}
} else {
if ($up_score >= $left_score) {
$matrix[$i][$j]{score}   = $up_score;
$matrix[$i][$j]{pointer} = "up";
}
else {
$matrix[$i][$j]{score}   = $left_score;
$matrix[$i][$j]{pointer} = "left";
}
}
}
}

##########################################################################################################
#                                                                       print matrix                                                                                                             #
########################################################################################################## 

my $j = length($seq1);
my $i = length($seq2);
my @seq1 = ("*");
my @seq2 = ("*");
my @s1=split(//,$seq1);
my @s2=split(//,$seq2);
push (@seq1,@s1);
push (@seq2,@s2);
if ($debug ==1 ) {

print "The Score Matrix\n";
print "-----------------------\n";
print "\t";
foreach (@seq1) {
print $_."\t";
}
print "\n";
for (my $u=0;$u<=$i;$u++) {
print $seq2[$u]."\t";
for (my $v=0;$v<=$j;$v++) {
print $matrix[$u][$v]{score}."\t";
}
print "\n";
}
 
#####################################################################################################
print "The Trace Matrix\n";
print "-----------------------\n";
print "\t";
foreach (@seq1) {
print $_."\t";
}
print "\n";
for (my $u=0;$u<=$i;$u++) {
print $seq2[$u]."\t";
for (my $v=0;$v<=$j;$v++) {
 
if ($matrix[$u][$v]{pointer} eq "none") {print  "N\t";}
if ($matrix[$u][$v]{pointer} eq "left") {print  "L\t";}
if ($matrix[$u][$v]{pointer} eq "up") {print  "U\t";}
if ($matrix[$u][$v]{pointer} eq "diagonal") {print  "D\t";}
 
}
print "\n";
}
 
}
my $align1 = "";
my $align2 = "";
 
my $score = $matrix[$i][$j]{score};
 
while (1) {
last if $matrix[$i][$j]{pointer} eq "none";
 
if ($matrix[$i][$j]{pointer} eq "diagonal") {
$align1 .= substr($seq1, $j-1, 1);
$align2 .= substr($seq2, $i-1, 1);
$i--; $j--;
}
elsif ($matrix[$i][$j]{pointer} eq "left") {
$align1 .= substr($seq1, $j-1, 1);
$align2 .= "-";
$j--;
}
elsif ($matrix[$i][$j]{pointer} eq "up") {
$align1 .= "-";
$align2 .= substr($seq2, $i-1, 1);
$i--;
}
}
 
$align1 = reverse $align1;
$align2 = reverse $align2;
 
print "Best Alignment\n";
print "--------------\n";
print "Seq1:\t$align1\n";
print "Seq2:\t$align2\n";
print "score = $score\n";
 
print time - $now." second\(s\)\n";

Last edited: 2011-12-26 23:04:06 +0100 (CET)

hlubenow

2011-12-26 23:57
User since
2009-02-22
875 Artikel
BenutzerIn
Hmm, sieht ja schonmal ganz gut aus.
Ich persönlich würde Code-Blöcke um 4 Leerzeichen einrücken, nach einer Weile der Umgewöhnung liest es sich dann besser (und ist daher inzwischen mehr oder weniger üblich geworden), z.B.
Code (perl): (dl )

1 2 3 4 5 6

for (my $i = 1; $i <= 10; $i++) { print "$i\n"; if ($i == 5) { print "i ist 5.\n"; } }
oder so.
Und sonst? Was ist die Frage?
topeg

2011-12-27 05:36

User since
2006-07-10
2611 Artikel
BenutzerIn

Der Code las sich als wäre er aus C portiert, darum habe ich ihn etwas geordnet und aufgeräumt.
Ich habe mal alle substr entfernt stattdessen zwei Arrays verwendet. Die Schleifen habe ich auch etwas verändert. Die Einrückungen habe ich nebenbei auch gemacht und die Möglichkeiten für Ein/Ausgabe über Dateien eingesetzt.
Insgesamt sollte der Code ein klein wenig schneller sein.
more (50.2kb)

View all threads created 2011-12-26 17:24.