Zufall in Matrix

[Eistoeter]

Hallo Leute,

ich programmiere schon sehr lange und sehr viel Java, aber gerade versuche ich ein Problem zu lösen, welches mich in den Wahnsinn treibt. Daher meine Hoffnung, dass jemand von euch vielleicht eine gute Idee hat.

Es geht um folgendes. Es gibt ein zwei-dimensionales NxN Array von double-Werten, das sieht z.B. so für N = 3 aus:

0.2 | 0.4 | 0.3
0.1 | 1.0 | 0.2
0.4 | 0.4 | 0.9

Ziel des Algorithmus ist es, aus jeder Reihe 1 Element auszuwählen. Die Zahlen geben dabei die Wahrscheinlichkeiten relativ zueinander an. Also eine Spalte mit 1.0 soll entsprechend wahrscheinlicher sein, als eine Spalte mit z.B. 0.4. Aber 1.0 heißt nicht, dass diese Spalte zu 100% gewählt wird. Es ist nur sehr viel wahrscheinlicher als eine Spalte mit 0.4.

Das ist eigentlich kein so schweres Problem. Ich habe erstmal so angefangen, dass ich die Summe jeder Reihe berechnet habe - und zwar bevor der eigentliche Algorithmus startet. Das wird also schon zu dem Zeitpunkt gemacht, bei dem das Array mit Werten gefüllt wird. Dann sah der Algo in etwa so aus:

Ansatz 1:

for (int i = 0; i < N; i++) {
   // Zufallszahl zwischen [0, Zeilensumme[
   double random = RANDOM.nextDouble() * rowSummations[i];

   // Über alle Felder gehen, dabei Summe mit hochzählen
   double sum = 0;
   for (int j = 0; j < N; j++) {
      sum += matrix[i][j];
      // Wenn gezählte Summe > Zufallszahl, dann Spalte gefunden
      if (sum > random) {
          result.setColForRow(i, j);
          break;
      }
   }
}

Nehmen wir z.B. die erste Zeile aus obigem Beispiel:
0.2 | 0.4 | 0.3 (Summe = 0.9)
-> Spalte 1 für Zufallszahlen [0 - 0.2[
-> Spalte 2 für Zufallszahlen [0.2 - 0.6[
-> Spalte 3 für Zufallszahlen [0.6 - 0.9[

So weit so gut. Jetzt ist dieser Algorithmus aber leider nicht besonders performant, im Worst Case O(n²). Mein Programm muss diesen Code aber mehr als 300 Millionen mal aufrufen. Und diese 300 Millionen mal sollten möglichst schnell gehen, da das Programm oft ausgeführt werden soll, um Auswirkungen verschiedener Parameter auf das Programmergebnis zu testen.
Mit anderen Worten -> dieses Problem muss irgendwie mit einem schnelleren Algorithmus gelöst werden.

N ist für die Daten mit denen ich arbeite praktisch immer 100. Eine 100 X 100 Matrix 300 Millionen mal ganz durchzulaufen ist schon viel. Also hatte ich folgende Idee:

Bevor der Algorithmus ausgeführt wird - zum Zeitpunkt der anfänglichen Belegung der Matrix mit den Werten - teile ich die Spalten in Gruppen auf. Bei N = 100 nehme ich 12 Gruppen zu jeweils 8 Spalten. 12 * 8 ist aber nur 96, deswegen ist die erste Gruppe 12 statt 8 Spalten groß. Zu jeder Gruppe speichere ich die Zwischensumme. So kann während des Algorithmus schnell bestimmt werden, zu welcher Gruppe die Zufallszahl gehört (Worst Case 12 Elemente durchgehen). Dann muss nur noch die Suche aus dem bisherigen Algorithmus pro Reihe auf maximal 8 Elemente (bzw. maximal 12 für die erste Gruppe) durchgeführt werden, anstatt auf 100. Das sollte schon einen großen Performance-Gewinn bringen. Der Algorithmus sah dann in etwa so aus:

Ansatz 2:

		for (int i = 0; i < N; i++) {
                        // Zufallszahl von [0 - Zeilensumme[  (die Gruppen-Summen sind in der zweiten Dimension von rowSummations)
			double randomNumber = RANDOM.nextDouble() * rowSummations[i][numberGroups - 1];

			// Finde richtige Gruppe
 			int group = 0;
			for (int j = 0; j < rowSummations[i].length; j++) {
				if (randomNumber < rowSummations[i][j]) {
					group = j;
					break;
				}
			}

			// Erster und letzter Eintrag von Gruppe bestimmen
			// Erste Gruppe ist potentiell größer als die anderen
			int first = group == 0 ? 0 : enlargementFirstGroup + colsPerGroup * group;
			int last = group == 0 ? colsPerGroup - 1 + enlargementFirstGroup : first + colsPerGroup - 1;

			// Suche Eintrag in Gruppe
			double sum = group == 0 ? 0 : rowSummations[i][group - 1];
			for (int j = first; j <= last; j++) {
				sum += matrix[i][j];
				if (randomNumber < sum) {
					result.setColForRow(i, j);
					break;
				}
			}
		}

Die zwei Algorithmen, die ich jetzt gepostet habe, funktionieren auch beide. Aber da ich es jetzt nicht 1-1 den alten Stand kopieren konnte, können sich kleine Fehler beim Tippen ins Forum ergeben haben. Darum geht es aber nicht wirklich.

Denn mein Problem kommt erst jetzt. Es gibt nämlich noch eine weitere Anforderung:
- Wurde eine Spalte in einer Reihe bereits ausgewählt, so darf sie in späteren Reihen nicht nochmals gewählt werden.

Dafür finde ich einfach keine Lösung - oder jedenfalls keine, die nicht die Performance-Verbesserung direkt wieder zunichte machen würde. Ich habe noch ein boolean[] ergänzt, das die Spalten repräsentiert und dort true ist, wo eine Spalte bereits ausgewählt wurde. Dann überprüfe ich das. Das große Problem dabei markiere ich im folgenden Code als Kommentar:

Ansatz 3:

                boolean[] addedCols = new boolean[N]; // NEU

		for (int i = 0; i < N; i++) {
                        // Zufallszahl von [0 - Zeilensumme[  (die Gruppen-Summen sind in der zweiten Dimension von rowSummations)
			double randomNumber = RANDOM.nextDouble() * rowSummations[i][numberGroups - 1];

			// Finde richtige Gruppe
 			int group = 0;
			for (int j = 0; j < rowSummations[i].length; j++) {
				if (randomNumber < rowSummations[i][j]) {
					group = j;
					break;
				}
			}

			// Erster und letzter Eintrag von Gruppe bestimmen
			// Erste Gruppe ist potentiell größer als die anderen
			int first = group == 0 ? 0 : enlargementFirstGroup + colsPerGroup * group;
			int last = group == 0 ? colsPerGroup - 1 + enlargementFirstGroup : first + colsPerGroup - 1;

			// Suche Eintrag in Gruppe
			double sum = group == 0 ? 0 : rowSummations[i][group - 1];
			for (int j = first; j <= last; j++) {
				sum += matrix[i][j];
				if (randomNumber < sum) {

                                        if (addedCols[j]) {

                                             // ---------------------------------------------
                                             // Bereits hinzugefügte Spalte erkannt, aber wie reagieren?
                                             //
                                             // Überlegung 1: So lange in Richtung des nächsten Nachbars gehen bis
                                             //    noch nicht verwendete Spalte erreicht
                                             //    -> Problem 1: Dann erhöht sich implizit die Wahrscheinlichkeit, dass die
                                             //                       benachbarten Spalten genommen werden relativ zu den
                                             //                       anderen noch zur Verfügung stehenden Spalten, obwohl an
                                             //                       den Wahrscheinlichkeiten nichts verändert werden soll
                                             //     -> Problem 2: Umso weiter am Ende des Algos wir sind, umso
                                             //                          unwahrscheinlicher wird es, direkt eine passende Spalte
                                             //                          zu finden (-> Performance-Gewinn wird wieder zu nichte gemacht)
                                             //
                                             // Überlegung 2: Zufallszahl-Reichweite einschränken, so dass die Werte bereits
                                             //                hinzugefügter Spalten von der Reichweite abgezogen werden.
                                             //     -> Problem 1: Wie performant realisieren - ohne für alle noch ausstehenden
                                             //                 Reihen mit einer Schleife die Subtraktions-Werte festzuhalten?
                                             //     -> Problem 2: Die Reichweiten der Gruppen stimmen dann nicht mehr. Da die
                                             //                 Zufallszahlen immer kleiner würden, würde man immer häufiger in den
                                             //                 niedrigen Gruppen landen (die Gruppen-Summen müssten also auch angepasst
                                             //                 werden -> Performance?)
                                             // ---------------------------------------------

				        }

					result.setColForRow(i, j);
                                        addedCols[j] = true; // NEU
					break;
				}
			}
		}

Wäre echt dankbar, falls jemand eine Idee hat?

Gruß

 

[kay73]

Ich hoffe mal, ich habe die Aufgabenstellung nicht total falsch verstanden. Abgesehen von Threads beim Iterieren der Zeilen sehe ich die einzige Optimiermöglichkeit auch nur bei der Spaltenwahl.

Nehmen wir z.B. die erste Zeile aus obigem Beispiel:
0.2 | 0.4 | 0.3 (Summe = 0.9)
-> Spalte 1 für Zufallszahlen [0 - 0.2[
-> Spalte 2 für Zufallszahlen [0.2 - 0.6[
-> Spalte 3 für Zufallszahlen [0.6 - 0.9[

Die etwas genaueren rechten Intervallgrenzen wären hier

[0-0.222(|(0.222-0,667)|(0,667-1,0]

, oder?

Nimm die als Zeilen einer Hilfsmatrix und mache darauf eine simple binäre Suche. Da kannst Du maximal O(log n) herausholen, indem Du zur Ur-Matrix eine Hilfsmatrix konstruierst, deren Zeilen eben die rechten Intervallgrenzen enthalten. Dann brauchst Du bei der Suche für N=100 im Schnitt 6-7 Tests im Gegensatz zu 50 bei linearer Suche.

Das Problem mit dem Verbot der Auswahl bereits ausgewählter Spalten ist mir nicht ganz klar: Nehmen wir an, wir sind in der dritten Zeile einer 10x10 Matrix und die Spalten 2 und 5 sind schon "verboten": Sollen dann die Spalten 1,3-4 und 6-10 erlaubt sein und mit Wahrscheinlichkeit relativ zu Ihrer Breite gewählt werden?

Einfacher: Angenommen im Beispiel oben sei in der ersten Zeile die Spalte 2 ausgewählt worden, wären die Spalten 1 (0.2) und 3 (0,3) noch erlaubt? Dann würde Spalte 1 mit Wahrscheinlichkeit 0,4 und Spalte 3 mit Wahrscheinlichkeit 0,6 gewählt?

 

[SlaterB]

wenn die Wahrscheinlichkeiten so schön in 0.1er Schritten bis maximal 1.0 oder auch 2.0 laufen, kann man gar alles auf int rechnen, 1-10 bzw. 1-20
und könnte komplette Arrays befüllen:
zu (0.2 | 0.4 | 0.3)
ein int[] x mit 1, 1, 2, 2, 2, 2, 3, 3, 3
und eine Zufallszahl z von 0-9 bzw. 0-8, dann ist x[z] direkt die errechnete Spalte,

zu viel sollte man davon nicht erwarten, nicht dass das ganze 10fach oder gar N-fach schneller wird,
allein die Bestimmung einer Zufallszahl könnte länger dauern als ein Array zu durchlaufen


mit Verbot von Spalten ist das aber leider alles hinfällig, nicht mal mehr die Summe kann gut im Voraus berechnet werden

 

[Eistoeter]

Sorry, das ich jetzt erst antworte (Ostern kam dazwischen). Vielen Dank für eure Anregungen.

Also erstmal: Die Wahrscheinlichkeiten sind eigentlich immer viel komplexer als so gerade, wie ich das jetzt in dem Beispiel gemacht habe. Da können die krummsten Kommazahlen rauskommen.

Das mit dem Verbot der Auswahl ist folgendermaßen: Das übergeordnete Ziel ist es, eine Zahlenfolge zu konstruieren. Also z.b. für N = 10 eine Zahlenfolge von 10 Zahlen (1-10). In dieser Folge darf jede Zahl nur einmal vorkommen. Die gewählte Spalte entspricht der Zahl, die verwendet wird. Also im Beispiel markiere ich mal die beispielsweise gewählten Spalten fettgedruckt:

0.2 | 0.4 | 0.3
0.1 | 1.0 | 0.2
0.4 | 0.4 | 0.9

Das führt zur Zahlenfolge: 2 1 3
Eine Zahlenfolge 3 3 1 wäre z.B. ungültig, weil die 3 zweimal vorkommt, dafür aber die 2 fehlt. Deswegen darf eine bereits ausgewählte Spalte nicht noch einmal ausgesucht werden. Eine verbotene Spalte sollte dann vom Algorithmus so behandelt werden, als würde sie gar nicht existieren.

@kay73: Du hattest ja eine Hilfsmatrix mit den Intervallgrenzen vorgeschlagen. So habe ich das im Prinzip ja umgesetzt. Die große Schwierigkeit ist dabei folgende.

Nehmen wir an wir haben die Beispiel 3X3 Matrix. Für die erste Zeile läuft es so:

1) Zufallszahl zwischen 0 und 0.9 generieren
2) Finde mit Binärsuche richtige Hilfsmatrix-Zeile
3) Suche Element in Zeile der Hilfsmatrix

So weit ist das noch alles gar kein Problem und auch schneller als der original Ansatz.

Für die zweite Zeile sieht es jetzt aber so aus (wenn wie oben die 2. Spalte in der ersten Zeile ausgewählte wurde):

1) Zufallszahl zwischen 0 und 0.3 generieren (nicht 1.3, da ja die 2. Spalte gesperrt ist und quasi nicht mehr existiert für uns)
2) Finde mit Binärsuche richtige Hilfsmatrix-Zeile
-> Problem: Die Intervallgrenzen stimmen jetzt nicht mehr, da unsere Zufallszahl nicht mehr so groß werden kann wie ursprünglich geplant. Hätten wir die Reichweite der Zufallszahl nicht verkleinert, könnte es sein dass wieder die Zufallszahl generiert wird, die uns genau in die Hilfsmatrix-Zeile der 2. Spalte bringt, die wir ja aber schon gewählt haben. Die Wahrscheinlichkeit für eine solche ungültige Zufallszahl wird immer größer mit der Zeit (da immer mehr Spalten gesperrt werden).
-> Lösung: Intervallgrenzen aller nachfolgenden Zeilen anpassen

Diese Lösung habe ich implementiert mit dem Ergebnis, dass die Performance sich leicht verschlechtert hat, gegenüber dem "naiven" Ansatz 1.

Kann es sein, dass es vielleicht für dieses Problem keine effizientere Lösung mehr gibt? Multithreading hilft mir nicht weiter, da ich das Programm insgesamt schon so gethreadet habe, dass alle CPUs auf 100% fahren.

 

[andiv]

Willst du eigentlich eine Lösung bei der die Wahrscheinlichkeiten aller gewählten Zahlen möglichst groß sind oder eine bei der die Zeilen stur nacheinander betrachtet werden und in der letzten Zeile dann halt der Wert genommen wird der noch übrig ist auch wenn er die kleinste Wahrscheinlichkeit hat? (Ich bin mir bewusst dass du nicht immer den größten Wert nimmst, sondern es nur am wahrscheinlichsten ist dass er genommen wird)

 

[Marco13]

Ist ja vielleicht eine ganz interessante Aufgabe. Aber nochmal zur Klärung:
- Die Matrix ist in der Größenordnung von ~100x100 (nicht 10000x10000 oder so...)
- Die Matrix wird EINmal gefüllt und dann nicht mehr geändert
- Auf Basis der Matrix werden dann SEHR oft neue "Indexverteilungen" berechnet

- D.h. die Aufgabe könnte von einem "RandomMatrixComputer" interface erfüllt werden, grob und NUR zur Verdeutlichung der Anwendung hingehackt sowas wie

import java.util.*;

interface RandomMatrixComputer
{
    void init(int rows, int cols, float matrix[][]);
    void compute(int result[]);
}

class RandomMatrixTest
{
    public static void main(String args[])
    {

        RandomMatrixComputer rmc = new DefaultRandomMatrixComputer();

        int rows = 50;
        int cols = 50;
        RandomMatrixTest r = new RandomMatrixTest(rows, cols, rmc);

        testSpeed(r, 10);
    }

    private static void testSpeed(RandomMatrixTest r, int runs)
    {
        long before = System.nanoTime();
        for (int i=0; i<runs; i++)
        {
            r.computeNewResult();
            System.out.println("result "+Arrays.toString(r.getResult()));
        }
        long after = System.nanoTime();
        System.out.println("duration "+(after-before)*1e-9);

    }

    private RandomMatrixComputer rmc;
    private int result[];

    public RandomMatrixTest(int rows, int cols, RandomMatrixComputer rmc)
    {
        this.rmc = rmc;
        Random random = new Random(1);
        float matrix[][] = new float[rows][cols];
        for (int r=0; r<rows; r++)
        {
            for (int c=0; c<cols; c++)
            {
                matrix[r][c] = random.nextFloat();
            }
        }
        rmc.init(rows, cols, matrix);

        result = new int[rows];
        Arrays.fill(result, -1);
    }

    public void computeNewResult()
    {
        rmc.compute(result);
    }

    public int[] getResult()
    {
        return result;
    }
}

Und jetzt geht es darum, dort einen "DefaultRandomMatrixComputer" zu erstellen, der die Aufgabe so schnell wie möglich löst?

Falls das alles korrekt ist, gibt es da sicher einige Ansätze :reflect: Welcher davon der "beste"/schnellste ist, läßt sich so im Voraus wohl nur schwer sagen.

Trotzdem ist die Vorgehensweise rein formal (!!!!) (zumindest für mich) im Moment noch fragwürdig. Ich weiß, dass man sich bei allem, was mit "Wahrscheinlichkeiten" zu tun hat, gewaltig auf die Fresse legen kann, wenn man "nach seinem Bauchgefühl" geht :autsch: (Ziegenproblem, Geburtstagsparadoxon, Simpsons-Paradoxon...). Deswegen ist alles, was ich sage, unter Vorbehalt zu sehen, weil ich kein Stochastik-Experte bin und das nicht alles im Kopf formal nachprüfen kann bzw. nicht mit Bleistift und Papier nachvollziehen will. Aber teilweise habe ich auch Bedenken wegen des Punktes, den andiv angedeutet hat: Die Wahrscheinlichheitsverteilungen der oberen Zeilen beeinflussen die Wahrscheinlichkeiten der unteren Zeilen - und die Verteilung in der letzten Zeile hat eine Aussagekraft von 0, weil einfach die letzte Spalte gewählt wird. Als suggestives Extrembeispiel, um den Punkt zu verdeutlichen:

1 0 0 0 0
0 1 0 0 0
0 0 1 0 0
0 0 0 1 0
1 1 1 1 0

In der letzten Zeile wird eine Spalte gewählt, die die Wahrscheinlichkeit 0 hat (oder nicht, oder es gibt eine Endlosschleife, oder eine Exception )

Je nachdem, worum es dabei geht, solltest du nochmal überlegen, ob das, was du da Programmieren willst (nochmal: rein formal) die Lösung des Problems ist. Es könnte sein, dass man die Wahrscheinlichkeiten alle als ganzes berücksichtigen und schon irgendwie miteinander verrechnen muss.

Falls ich das alles richtig sehe, wäre eine formal richtige Verteilung in diesem Sinne eine, die so etwas macht wie "Collections.shuffle" auf dem Array der Spaltenindizes, wobei aber die Positionen der indizes am Ende nicht normalverteilt sind (wie bei shuffle) sondern irgendwie von irgendwelchen relativen oder bedingten Wahrscheinlichkeiten abhängen...

 

[Eistoeter]

Ich will eine Lösung, in der jede Zeile stur nacheinander betrachtet wird. Es stimmt, mit jeder Auswahl in einer Zeile werden die Wahrscheinlichkeiten der unteren Zeilen beeinflusst, so dass die letzte Zeile eigentlich gar keine Relevanz mehr hat, sondern dort einfach nur noch das gewählt wird, was übrig bleibt - das ist so gewollt.

Einen Punkt hätte ich vielleicht doch noch sagen sollen: Die Matrix verändert sich schon nach jedem Durchlauf, sie wird nicht einmal initialisiert und bleibt dann immer gleich. Nach jedem Durchlauf werden die Wahrscheinlichkeiten aktualisiert (alle gewählten Spalten erhalten einen Bonus wenn die Lösung gut war, wohingegen alle anderen einen Malus bekommen). Das ganze ist aber auf Minimal- und Maximal-Werte begrenzt (z.B. 0.1 und 1).

Wer's noch genauer wissen will: Es handelt sich um das Max-Min-Ant-System. Die Matrix stellt Pheromon-Werte da, die von Ameisen hinterlassen wurden. Bei jedem Durchlauf startet eine neue Ameise von oben und läuft durch die Matrix. Dort wo am meisten Pheromon ist, soll sie natürlich am ehesten hin gehen. Dabei besucht sie keinen Ort zweimal. Die Sache mit dem "in der letzten Zeile bleibt nur noch eine Auswahl übrig" ist somit OK, da die Ameise dann halt keine anderen Möglichkeiten mehr hat. Es geht nur darum, bei jedem Entscheidungspunkt von den möglichen Entscheidungen die Wahrscheinlichkeiten richtig auszuwerten.

 

[Marco13]

OK, wenn die Matrix sich bei jedem Schritt ändert macht das ja alles, was man sich da an Vorberechnungen oder geschickten Datenstrukturen ausdenken könnte wieder zunichte. Zumindest kann ich mir kaum vorstellen, wie man da dann strukturell oder direkt algorithmisch noch etwas rausholen könnte. Die Sache mit den schon verwendeten Spalten könnte man dann vielleicht noch etwas geschickter machen, z.B. die Spalten in eine Collection legen, und die schon benutzten einfach raussrteichen, aber VORberechnungen von Summen machen ja keinen Sinn, wenn sie sowieso nur einmal verwendet werden. Darüberhinaus könnte man überlegen, ob das nicht schön in CUDA oder OpenCL umzusetzen wäre, aber das ist ja nochmal eine ganz andere Frage.

 

[Eistoeter]

Also die Matrix ändert sich nicht nach jedem Schritt (Zeile), sondern nach jedem kompletten Durchlauf.

Die Summen zu berechnen ist mehr oder weniger kostenlos, da dies bei der Aktualisierung der Matrix geschehen kann, bei der sowieso die ganze Matrix durchgegangen werden muss.

Aber mir wird das ehrlich gesagt alles zu schwammig. Ich glaube ich lebe lieber mit dem ersten Ansatz (der ist wenigstens einfach zu verstehen) und suche an anderen Stellen nach Performance-Verbesserungen.

Danke auf jeden Fall für eure Überlegungen!

 

[Marco13]

Es ging darum, dass man sich "ausgefeiltere" Dinge hätte überlegen können, wenn diese dann nicht auch 300 Millionen mal hätten gemacht werden müssen.
BTW: Die Summen hängen doch davon ab, welche Spalten vorher besucht wurden, also können sie auch nicht direkt vorher berechnet werden...?
Aber egal, wenn du einen anderen Weg suchen willst...

 

[Eistoeter]

Ja genau, die Summen hängen davon ab, welche Spalten vorher besucht wurden. Zumindest wenn man es auf die Art löst, das man die Summen nach jeder Zeile reduziert um eine korrekte Zufallszahl zu erhalten. Das ist auch genau das Problem, warum der Ansatz mit der Hilfsmatrix nicht schneller ist.

Aber es kommt ja auch noch dazu, dass nach jedem Durchlauf ein Matrix-Update stattfindet, in dem alle Werte erhöht / erniedrigt werden. Dort können die initialen Summen berechnet werden.

Wenn man einen Weg finden würde, der nicht erfordert die Summen nach jeder Zeile zu reduzieren, das wäre was. Aber da renn ich vor ne Wand.