Seit dem letzten Post ist viel los, einschließlich der Vorbereitungen JMP 16 , das erst vor ein paar Tagen herauskam, und ohne Antworten auf meinen letzten Post, in dem ich nach Spaltenformel-Herausforderungen fragte, auf die Leser gestoßen waren, dachte ich, ich hätte nicht die Zeit oder das Material, um etwas zu posten, das es wert ist, gelesen zu werden.
Ich hatte in letzter Zeit etwas mehr Zeit und es ist Zeit für eine neue Herausforderung. Wir stellen die Spaltenformel-Herausforderungen ein anderes Mal auf und beginnen mit einer neuen Herausforderung und einem neuen Ansatz. In diesem Beitrag und wahrscheinlich in einigen danach werde ich darüber sprechen eine allgemeine Frage, die auf einer Antwort von Craige Hales auf einen Community-Beitrag basiert . Craige ( @Craige_Hales ), obwohl er sich von JMP Development zurückgezogen hat, bleibt hier in der Community sehr aktiv, beantwortet Benutzerfragen und pflegt seine Uncharted-Blog , und Beiträge zu diesem Blog. Die Frage war nach der Effizienz, eine (mehrdimensionale) Liste ohne Verwendung von Schleifen zu transponieren. Ich möchte eine allgemeinere Frage stellen, die ich mir oft gestellt habe: „Gibt es bessere Alternativen zur Iteration, und wenn ja, welche?“ In den meisten Fällen bedeutet besser schneller, aber ich werde auch die Kompaktheit des Codes berücksichtigen, wenn die Geschwindigkeit nicht zu sehr beeinträchtigt wird. Wir werden uns auch verschiedene Iterationsansätze ansehen, um zu sehen, ob es auch dort Unterschiede gibt. Beginnen wir mit der Textsuche.
Saiten
Eine häufige Aufgabe in der Bioinformatik besteht darin, eine Zeichenfolge zu durchsuchen, die genetische Informationen enthält, um bestimmte Muster von Teilzeichenfolgen zu finden. Bei DNA bestehen die Zeichenketten aus nur vier Buchstaben: A, T, C und G. Die Datei Covid19.TXT (an diesen Beitrag angehängt) enthält die ungefähr 30.000 Zeichen lange genetische Sequenz für das Covid19-Genom (im Gegensatz dazu ist das menschliche Chromosom 1, das längste, ungefähr 290 Millionen Zeichen lang). Wie schnell kann gesucht werden, um den Ort aller Vorkommen bestimmter Teilzeichenfolgen mit einer Länge von 2 bis 10 Zeichen zu finden?
Sehen wir uns vier Ansätze an, drei davon iterativ. Beachten Sie unabhängig von der Technik, dass Übereinstimmungen von Teilzeichenfolgen überlappen können. Zum Beispiel, BAB ist doppelt enthalten BABAB , von den Positionen 1 bis 3 und wieder von 3 bis 5. Wir brauchen etwas, das diese beiden Ereignisse erfasst. Wir beginnen mit Substr Verwenden Sie den folgenden Code, wo inString ist die zu durchsuchende Zeichenfolge und inSeq ist die Teilzeichenfolge, nach der gesucht werden soll.
findSequence = Function({inString,inSeq},{i,seqLen,checkTo,outVec},
seqLen = Length(inSeq);
checkTo = Length(covid19)-seqLen+1; //Only need to go until the substring includes the last character
outVec = J(checkTo,1,0);
For(i=1,i<=checkTo,i++,
outVec[i] = Substr(inString,i,seqLen) == inSeq;
);
Loc(outVec);
);
Dies war der am wenigsten beeindruckende Ansatz und dauerte durchschnittlich 1,2 Sekunden, unabhängig von der Länge der Teilzeichenfolge. Die Ergebnisse sind unten gezeigt (in Mikrosekunden). Der Test generierte 10 zufällige Teilstrings für jede Teilstringlänge und wiederholte die Suche 20 Mal. Die Werte wurden über alle Teilstrings einer gegebenen Länge zusammengefasst. Vollständige Details für einzelne Teilstrings sind in Covid 19 Ergebnisse Vollständig.JMP an diesen Beitrag angehängt. Der zum Generieren des Tests verwendete Code befindet sich in Test Covid19 Search.JSL.
| Teilstringlänge |
Substr |
Index |
Enthält |
PatMatch |
| 2 |
1.227.342 |
60.312 |
79.640 |
4.046 |
| 3 |
1.245.909 |
62.650 |
22.143 |
3.103 |
| 4 |
1.215.641 |
66.045 |
4.938 |
2.881 |
| 5 |
1.218.983 |
72.337 |
1.293 |
2.880 |
| 6 |
1.220.846 |
76.705 |
357 |
2.956 |
| 7 |
1.249.843 |
77.127 |
234 |
3.010 |
| 8 |
1.212.982 |
84.413 |
132 |
3.078 |
| 9 |
1.227.892 |
87.479 |
113 |
3.331 |
| 10 |
1.208.948 |
89.732 |
106 |
3.742 |
Das haben wir bei früheren Herausforderungen gesehen Substr ist langsam. Was wäre, wenn wir den String „atomisieren“, ihn in einzelne Zeichen zerlegen und Zeichenbereiche betrachten würden? Die Motivation dafür liegt in der Größenordnung des direkten Zugriffs auf Zeichenfolgenwerte über die Indizierung, wie Sie es in einer Sprache wie C können. Der Atomisierungsteil wird in JSL benötigt, da Sie eine Zeichenfolge nicht direkt indizieren können. Der Code, den ich verwendet habe, ist hier:
findSequence = Function({inString,inSeq},{i,seqLen,checkTo,outVec,inList,seqList},
seqLen = Length(inSeq);
checkTo = Length(inString)-Length(inSeq)+1; //Only need to go until the substring includes the last character
outVec = J(checkTo,1,0);
inList = Items(inString,"");
seqList = Items(inSeq,"");
For(i=1,i<=checkTo,i++,
outVec[i] = inList[i::(i+seqLen-1)] == seqList;
);
Loc(outVec);
);
Dieser Ansatz scheint zwar zeitaufwändiger zu sein, da zusätzlich zur Suche der anfängliche Aufwand für das Aufbrechen der gesuchten Zeichenfolge anfällt, er ist jedoch etwa 13-19-mal schneller, wie in der Ergebnistabelle zu sehen ist. nicht wie Substr , jedoch erhöht sich die Suchzeit als Funktion der Teilstringlänge.
Der letzte iterative Ansatz verwendet Enthält . Die Idee besteht darin, die angegebene Teilzeichenfolge zu finden, ihre Position in der gesuchten Zeichenfolge zu markieren und den Vorgang ein Zeichen nach der Position zu wiederholen, bis die Teilzeichenfolge nicht gefunden werden kann. Der Reiz dieser Methode besteht darin, dass die Suche sehr schnell wird, wenn die Teilzeichenfolge in der gesuchten Zeichenfolge nur wenige Vorkommen hat. Der folgende Code ist Brady zu verdanken Brady ( @brady_brady 
findSequence = Function({inString,inSeq},{pos=0,ans={}},
while(1,
pos = contains(inString, inSeq, pos);
if(pos, insertinto(ans, pos++), break())
);
ans;
);
Wie Sie sehen können, beginnt dieser Ansatz zwar etwas langsamer als die Indizierungsmethode für zweistellige Teilzeichenfolgen, wird aber schnell schneller, wenn die Teilzeichenfolgen länger werden, dh, da weniger Teilzeichenfolgen zu finden sind.
Der letzte Ansatz verwendet das, was zu einem meiner Lieblings-JSL-Tools wird, den Musterabgleich. Der Code ist nicht nur schnell, sondern auch kurz.
findSequence = Function({inString,inSeq},{i=1,locs={}},
Pat Match(inString,Pat Pos()>>r + inSeq + Pat Test(Insert Into(locs,r+1);0));
locs;
);
Während es für die zweistelligen Teilzeichenfolgen etwas langsamer ist, scheint es nicht von der Länge der Teilzeichenfolge für drei Zeichen und mehr beeinflusst zu werden.
Ganz klar, der iterative Ansatz mit Enthält und der Mustervergleichsansatz sind überlegene Methoden zum Suchen langer Zeichenfolgen. Aber welche Methode und wann? Da die Tests mit verschiedenen Teilzeichenfolgen durchgeführt wurden, gibt es Hinweise darauf, wie sich die Häufigkeit der Teilzeichenfolgen (Anzahl Vorkommen*1000/29.903) auf die Suchzeit auswirken kann. Das folgende Diagramm zeigt einen Kreuzungspunkt irgendwo um eine Häufigkeit von 2 Teilstrings pro 1.000 Zeichen herum. Da alle Tests mit demselben covid19-Genomstring durchgeführt wurden, ist jedoch unklar, wie die Suchzeit von der Länge des gesuchten Strings beeinflusst wird.

Um die Antworten auf diese Fragen zu finden und mögliche Wechselwirkungen zwischen Faktoren zu untersuchen, habe ich ein Experiment mit einem raumfüllenden Design durchgeführt, bei dem die Länge der gesuchten Zeichenfolge (10.000–500.000 Zeichen), die Länge der Teilzeichenfolge (2–21 Zeichen) und die Häufigkeit ( 0,75 – 3 Teilstrings pro 1.000 Zeichen). 125 Läufe wurden unter Verwendung eines Fast-Flexible-Designs mit gerundeten Bruchteil-String-Längenwerten generiert. Für jeden Durchlauf wurde eine zufällige Folge alphabetischer Zeichen generiert. Innerhalb dieser Zeichenfolge wurde an zufälligen Stellen eine Teilzeichenfolge eingefügt, die aus zufällig generierten Zahlen bestand. Da Häufigkeitswerte zu einer nicht ganzzahligen Anzahl von Einfügungen führen könnten, wurde die Anzahl der Einfügungen als Obergrenze der Häufigkeit multipliziert mit der Länge der gesuchten Zeichenfolge geteilt durch 1.000 genommen. Jeder Durchlauf wurde 10 Mal mit der gleichen Kombination aus gesuchter Zeichenkette/Teilzeichenkette wiederholt. Der Code, der zum Generieren der Läufe verwendet wird, ist in Führen Sie Design Script.JSL aus und Ergebnisse sind da Lückenfüllungsergebnisse 1.JMP an diesen Beitrag angehängt.
Eine schnelle Überprüfung der Ergebnisse zeigt nur 13 Beobachtungen, bei denen die Enthält Methode übertrifft den Musterabgleich. Von diesen ist der größte Unterschied etwa 4.000 Mikrosekunden schneller (während der Musterabgleich ganze 0,72 Sekunden schneller war als Enthält ). Das Anpassen eines Gaußschen Prozessmodells an die Daten legt nahe, dass die Länge der gesuchten Zeichenfolge den größten Effekt mit einer starken Häufigkeit durch die Wechselwirkung der Länge der gesuchten Zeichenfolge hatte. Die Länge der Teilzeichenfolge hatte wenig Einfluss. Anfänglich wurde ein lineares Modell an die Daten angepasst, aber die Residuenanalyse legte nahe, dass die Anpassung fragwürdig sein könnte. Ich habe gepostet eine interaktive Version des Profilers für dieses und die Nachfolgemodelle an JMP Public.
Ein zweites raumfüllendes Design wurde erstellt, um eine bessere Darstellung von Fällen zu erhalten, in denen Enthält Musterabgleich übertroffen. Für dieses Experiment blieben die Bereiche für beide Saitenlängen wie zuvor und die Frequenz wurde zwischen 0 und 1 variiert. Die Ergebnisse sind da Lückenfüllende Ergebnisse 2.JMP . Die Richtung der Ergebnisse war ähnlich wie zuvor, wobei längere gesuchte Zeichenfolgenlängen und höhere Frequenzen den Musterabgleich begünstigen.
Das Kombinieren beider Ergebnissätze und das Anpassen eines Gaußschen Prozessmodells legt nahe, dass für Zeichenfolgen von 50.000 Zeichen Enthält ist überlegen, wenn die Frequenzen unter 2 liegen. Bei 100.000 Zeichen sinkt die Häufigkeit auf 1 und bei 200.000 Zeichen unter 0,5. Die Schlussfolgerung daraus ist, dass Sie für spärlichere kürzere Saiten verwenden Enthält und wenn die Stringlänge und -dichte zunimmt, wird der Musterabgleich überlegen sein.
Herausforderung des Monats
Finden Sie mithilfe des Covid19-Genoms die längste Teilzeichenfolge, die zwei- oder mehrmals vorkommt.
Anstatt dies bis zum nächsten Beitrag unbeantwortet zu lassen, lassen Sie uns diskutieren und teilen. Ich glaube, das ist ein sehr schwieriges Problem, und wir werden von einer offenen Diskussion profitieren.
Inklusive Überschneidungen? Das macht es einfach. Ich bin mir nicht sicher, wie ich es machen würde, ohne die Überlappungen. Diese Lösung ist eine 5-10-Sekunden-Antwort für die gelieferten Daten. Und es würde für viel längere Daten mit viel längeren Teilzeichenfolgen schmerzhaft werden.
Die Idee ist, zu zählen, wie viele von jedem Teilstring für verschiedene Teilstringlängen an jeder Position 1 ... Textlänge vorkommen. Fragen Sie dann "Gibt es Zählungen > 1?" Verwenden Sie eine binäre Suche nach der Teilstringlänge, um die größte Teilstringlänge zu finden, die Duplikate erzeugt. Am Ende sollte ein wenig mehr Arbeit geleistet werden, um zu sehen, ob es einen anderen, interessanteren Teilstring derselben Länge gibt.
txt = Load Text File( "z:/covid-19.txt" );
// finding the longest duplicated substring, including overlaps.
// if there is one, then there will be shorter ones of every possible
// shorter length. So, as the possible length of the substring grows,
// there is a boundary at some length where all shorter lengths have
// duplicates and no longer lengths have duplicates. find the boundary.
// do a binary search for the longest string that occurs 2 or more times
// testlength is the length during a pass
// see https://en.wikipedia.org/wiki/Binary_search_algorithm#Procedure
// for a *real* binary search. this code is potentially wrong and very
// poorly tested!
// Changing the maxlength guess below from 1...1000 will produce different times.
// 20 is a good starting point, bigger than 200 or so is very slow because
// length(txt) * maxlength strings are being manipulated, a lot.
start = Tick Seconds();
minlength = 1;
maxlength = 27; // initial guess, will (slowly) auto-adjust bigger as needed
firstpass = 1; // needed to detect guess too small
minkey = .;
minval = .;
while( 1,
testlength = if( firstpass == 1,
Min( Length( txt ), maxlength )// this binary search starts at the big end, then splits on subsequent passes. This lets code
, // else
Floor( (minlength + maxlength) / 2 )// below test that the big end is big enough on the first pass.
);
Write( "testing strings of length ", testlength, "\!n" );
Wait( 0 );
aa = [=> 0];
patmatch( txt,
Pat Len( testlength ) >> sequence + //
pattest( // build the associative array of key=sequence, val=count
aa[sequence]++;
// usually I use 1; on the next line. that always succeeds.
// here, I want pattest to make the matcher back up and retry...
0; // fail! try again at next position
),
NULL, // unused positional replacement
FULLSCAN // needed to bypass a heuristic and get the last little bit
);
keys = aa << getkeys;
vals = aa << getvalues;
idx = Loc Max( Matrix( vals ) );
bigN = vals[idx];
//
// the binary search is for the biggest testlength that has bigN > 1.
// bigN is the maximum repeats of a string of testlength characters.
//
if( bigN == 1,
maxlength = testlength;// maxlength promises to be 1-too-big at the end
, //else
minlength = testlength;
Write( " found some, minlength=", minlength, " maxlength=", maxlength, "\!n" );
minkey = keys[idx];
minval = vals[idx];
if( firstpass == 1, // looking for a big-enough high end, still
maxlength = Min( 2 * maxlength, Length( txt ) );
Write( "Increasing max, minlength=", minlength, " maxlength=", maxlength, "\!n" );
Continue();
);
);
firstpass = 0;
if( minlength >= maxlength - 1,
if( Is Missing( minkey ),
Continue();// needs one more pass
, // else
Break();// got it
)
);
);
Write( "\!n", minkey, " occurs ", minval, " times at " );
patmatch(
txt,
Pat Pos() >> location + //
minkey + //
pattest(
Write( location + 1, " " );
1;
) + //
Pat Fail() // this will also retry the match, 1 character further along
);
Write( "\!n" );
stop = Tick Seconds();
Show( stop - start );