Génération de matrices sous contraintes

En python, sur chaque ligne et chaque colonne

a marqué ce sujet comme résolu.

Bonjour les zestotifiques ! :)

Pour un petit problème que j’ai imaginé, je dois générer des matrices qui respectent des contraintes sur chaque ligne et chaque colonne. J’avais une idée d’algorithme qui me semblait marcher mais impossible de le mettre en oeuvre en python.

Je viens donc vous quémander de l’aide, car sinon je vais rester bloqué et je n’ai pas envie d’écrire ces matrices à la main. :(

Mon problème donc, je dispose de mes contraintes (toujours des nombres entiers strictement positifs) :

  • un tableau de contraintes pour les lignes : [l0,l1...lm1l_0,l_1... l_{m-1}]
  • un tableau de contraintes pour les colonnes : [c0,c1...cn1c_0,c_1...c_{n-1}]
  • avec m>=nm>=n

Et j’aimerais générer une matrice telle que:

  • Pour chaque ligne ii,la somme de la ligne ii est égale à lil_i
  • Pour chaque colonne ii, la somme de la colonne ii est égale à cic_i

Ah oui, il faudrait que les éléments de ma matrice semblent aléatoires, donc pas juste une suite de 0 puis lil_i par exemple.

Pour ça j’avais pensé à l’algorithme suivant :

Pour i=0 à m-1:
   l' <- séquence aléatoire de m-1-i éléments avec leur somme = l_i
   matrice[i][i:] = l'
   contrainte_colonne = contrainte_colonne - l'
   c' <- séquence aléatoire de n-2-i éléments avec leur somme = c_i
   matrice[:][i+1:] <- c'
   contrainte_ligne = contrainte_ligne - c'

Malheureusement mon implémentation ne fonctionne pas à cause de mes séquences qui ne font pas la même taille entre autres.

Voici un extrait en python 3 avec numpy :

    for i in range(m):
        line = np.array(random_ints_with_sum(contrainte_ligne[i],max_length=m-1-i))
        mat[i][i:] = line
        contrainte_ligne -= line
        col = np.array(random_ints_with_sum(contrainte_colonne[i],max_length=n-2-i))
        mat = mat.transpose()
        mat[i][i+1:] = col
        mat = mat.transpose()
        contrainte_colonne -= col

Est ce que vous auriez des idées pour améliorer cet algorithme / le corriger / en trouver un qui marche mieux ?

Merci d’avance !

+0 -0

Salut,

À ta place je créerais une matrice (m1)×(n1)(m-1)\times (n-1), remplie de nombres entiers aléatoires positif, puis je complèterai la colonne et la ligne manquantes avec des nombres (possiblement negatifs) qui respectent les contraintes. Les deux dernières contraintes donnent le terme dans le coin final (si une grille vérifiant les contraintes existe, ce qui n’est pas le cas général d’ailleurs, il faut déjà que la somme des contraintes soit la même dans les deux directions). Ensuite il ne reste plus qu’à redistribuer les valeurs négatives en repércutant les changements à la fois sur les lignes et les colonnes pour ne pas perdre le respect des contraintes. Si c’est impossible à un moment de rester positif partout, c’est que la grille n’a pas de solution…

+1 -0

Peu importe pour le fais que ça semble aléatoire ou non. Au pire, à la fin tu peux mélanger les lignes et les colonnes.

Après, pour une suite de 00 puis lil_i oui bon ça aide pas de mélanger.

+0 -0

Au pire, à la fin tu peux mélanger les lignes et les colonnes.

Heu… Comment tu fais ça ? A priori les lil_i et cic_i sont différents les uns des autres…

+0 -0

Échanger deux colonnes ne devrait pas poser de problèmes si ?
Pareil pour les lignes ?

En échangeant aléatoirement des lignes on les mélanges et pareils ensuite pour les colonnes ?

Ah oui, je vois, ça échange également les contraintes des deux lignes échangées. (╥_╥) Pareil pour les colonnes. :/

+0 -0

Il y a des contraintes sur les contraintes :

  • la somme des Li est nécessairement égale à la somme des Ci.
  • chaque Li doit être supérieur ou égal à n
  • chaque Ci doit être supérieur ou égal à m

Il faut commencer par vérifier que les contraintes peuvent être respectées.

Ensuite, on peut partir d’une matrice triviale remplie de 1.
on calcule Sc le tableau des sommes des colonnes et Sl le tableau des sommes des lignes. On cherche le plus petit élément de Sc qui vérifie la contrainte, pareil pour Sl. On ajoute 1’élément correspondant, et on recommence.

On s’arrête quand toutes les contraintes sont respectées.

+0 -0

Si pour chaque élément de la matrice on fait v[i, j] = l[i]*c[j]/total on arrive à une solution valide mais avec des nombres flottants, après il faut gérer l’arrondi.

Exemple avec 5 lignes / 4 colonnes

Sommes de chaque ligne
30
22
24
19
27

Sommes de chaque colonne
36 26 31 29

Total 122

Valeurs calculées
8.852 6.393 7.623 7.131
6.492 4.689 5.590 5.230
7.082 5.115 6.098 5.705
5.607 4.049 4.828 4.516
7.967 5.754 6.861 6.418

EDIT: bon ça donne des valeurs plutôt "pondérées", si on veut que cela paraissent plus aléatoire on peut éventuellement appliquer un offset +/- à chaque valeur dont la somme est nulle en ligne / colonne.

+1 -1

Normalement en utilisant le formalisme de la transformée Mojette.

Les contraintes deviennent qu’à partir des projections selon le vecteur (0,1)(0, 1) (lnl_n) et le vecteur (1,0)(1, 0) (cmc_m).

Ils nous manque bien sûr des projections (sauf dans le cas où n = m = 2) puisque selon le lemme de Katz la convolution des fantômes correspondant tiens dans le support (la matrice).

Je me demande si ajouter des projections jusqu’à que la convolution des fantômes des vecteurs remplisse l’espace du support donnerait un résultat cohérent.

Je suppose que oui mais il faudrait avoir certaines conditions pour que le résultat soit positif, sans quoi et bien il ne le sera pas :/

+0 -0

Et bien et bien, je m’attendais pas à un tel enthousiasme !

Je ne vais pas pouvoir répondre à tout les messages par contre je suis désolé :/

En regardant à tête reposée le problème je me suis aussi dit que peut être c’était simplement un problème de programmation linéaire, donc peut être résoluble avec la méthode du simplexe, mais c’est sûrement de l’artillerie lourde pour pas grand chose.

Pour ce qui est de mélanger les lignes et les colonnes, ça devrait passer en faisant bouger des lignes et colonnes entières avec leurs contraintes associées : Si on intervertit la ligne i avec la ligne j on intervertit aussi les contraintes lil_i et ljl_j dans le tableau et idem pour les colonnes.

Il y a des contraintes sur les contraintes :

  • la somme des Li est nécessairement égale à la somme des Ci.
  • chaque Li doit être supérieur ou égal à n
  • chaque Ci doit être supérieur ou égal à m
etherpin

Effectivement je n’avais pas pensé à ces contraintes (quand on y pense c’est vrai qu’elles sont immédiates pourtant ^^' ), ma génération de contraintes ne les prenait pas en compte et pourtant ça marchait, sûrement un coup de chance

Ensuite, on peut partir d’une matrice triviale remplie de 1.
on calcule Sc le tableau des sommes des colonnes et Sl le tableau des sommes des lignes. On cherche le plus petit élément de Sc qui vérifie la contrainte, pareil pour Sl. On ajoute 1’élément correspondant, et on recommence.

On s’arrête quand toutes les contraintes sont respectées.

etherpin

Alors là j’avoue que je ne comprend pas l’algorithme, autant celui de @ache je vois l’idée même si la partie redistribution me semble plus compliquée à mettre en place.

Si pour chaque élément de la matrice on fait v[i, j] = l[i]*c[j]/total on arrive à une solution valide mais avec des nombres flottants, après il faut gérer l’arrondi.

fromvega

Le problème va justement être avec l’arrondi je pense, on n’est obligé de gérer les deux contraintes pour chaque case en parallèle. Du coup à chaque fois les valeurs risquent de potentiellement se propager dans toute la matrice et ça devient tout de suite un algorithme avec une forte complexité.

Normalement en utilisant le formalisme de la transformée Mojette.

Les contraintes deviennent qu’à partir des projections selon le vecteur (0,1)(0, 1) (lnl_n) et le vecteur (1,0)(1, 0) (cmc_m).

Ils nous manque bien sûr des projections (sauf dans le cas où n = m = 2) puisque selon le lemme de Katz la convolution des fantômes correspondant tiens dans le support (la matrice).

Je me demande si ajouter des projections jusqu’à que la convolution des fantômes des vecteurs remplisse l’espace du support donnerait un résultat cohérent.

Je suppose que oui mais il faudrait avoir certaines conditions pour que le résultat soit positif, sans quoi et bien il ne le sera pas :/

ache

Là je suis complètement perdu mais je vais regarder un peu ce dont tu parles ça me semble intéressant !

Finalement j’ai trouvé une solution, un peu idiote, mais qui fonctionne en changeant mon problème.

Dans mon problème je devais générer des instances composées des triplets suivants :

  • Contrainte lignes
  • Contrainte colonnes
  • Matrice sous contrainte

Le problème se posait du fait que je générais les contraintes d’abord, puis tentait de générer la matrice ensuite. Mais si je génère juste les contraintes de lignes, puis la matrice qui les satisfait (chaque ligne aléatoirement qui respectent la contrainte ligne) alors les contraintes colonnes sont elles aussi générées automatiquement (il suffit de sommer les colonnes de la matrice).

Finalement c’était juste une question de point de vue. ;)

Par contre je vais continuer à regarder le problème de la génération avec tout ce que vous m’avez dit, ça me semble intéressant à poser comme autre problème !

+0 -0

Pour ce genre de problème (générer aléatoirement uniformément une solutions à un problème), une stratégie qui marche bien en théorie et en pratique est de définir des opérations de perturbation locales qui transforment une solution valide en une autre solution. L’algorithme est alors le suivant : on part d’une solution à notre problème (pas du tout aléatoire), on tire aléatoirement une "petite" perturbation, et on répète ça un très grand nombre de fois.

Autrement dit, on définit un graphe sur l’ensemble des solutions, on part d’un nœud quelconque, et à chaque étape on choisit aléatoirement un de ses voisins. En fait, on peut montrer que si ce graphe satisfait certaines conditions (connexité, symétrie, …), alors, quand le nombre d’étapes tend vers l’infini, on a autant de chance d’être sur chaque sommet : la loi est uniforme sur l’ensemble des solutions (plus exactement, converge vers la loi uniforme, en un certain sens).

Si tu as déjà vu une chaîne de Markov (qui n’est rien d’autre qu’un graphe avec des probabilités écrites sur les arêtes), tout ceci n’est pas un miracle : une chaîne de Markov réversible, ergodique, symétrique a pour unique loi stationnaire la loi uniforme. Avec ce formalisme, on est intéressé par des versions quantitatives de la convergence décrite plus haut — le temps de mélange, à savoir le plus petit temps tt tel que pour n’importe quel état de départ, à partir de l’étape tt, la distance entre notre distribution sur l’ensemble des sommets et la loi uniforme est inférieure à 0.01 (au sens de la distance en variation totale, par exemple). C’est une manière de fixer le nombre d’étapes à effectuer en fonction de à quel point tu t’autorises à dévier de la loi uniforme sur l’ensemble des solutions.

Evidemment, toute la difficulté du problème consiste à définir ces opérations de perturbation, de manière à ce qu’on puisse aller de n’importe quelle solution valide à n’importe quelle autre par une suite de perturbations, et que la chaîne de Markov converge "rapidement". Pour ce problème de matrices sous contrainte, connu — je découvre également — sous le nom de "contingency tables" dans la littérature, il y a une chaîne de Markov classique, apparemment due à Diaconis et Gangolli :


Pour perturber la matrice MM :

  • Générer uniformément aléatoirement 1i,im,1j,jn1 \le i, i'\le m, 1\le j, j' \le n tels que iii \neq i' et jjj\neq j'.
  • Modifier MM de la manière suivante : Mij=Mij+1M_{ij}=M_{ij}+1, Mij=Mij1M_{i'j}=M_{i'j}-1, Mij=Mij1M_{ij'}=M_{ij'}-1, Mij=Mij+1M_{i'j'}=M_{i'j'}+1.
  • Si une des cases de MM devient négative, revenir à l’ancien MM (ne rien changer).

Pour générer une matrice sous contrainte (quasi-uniformément) aléatoirement :

  • Partir de n’importe quelle matrice MM respectant les contraintes
  • Répéter un très grand nombre de fois :
    • Perturber MM

Pour choisir ce "très grand nombre de fois", j’imagine que tu peux le faire expérimentalement. Les meilleurs résultats théoriques sur le temps de mélange de la chaîne de Markov de Diaconis et Gangolli sont à mon avis assez loin de la réalité.

Si pour chaque élément de la matrice on fait v[i, j] = l[i]*c[j]/total on arrive à une solution valide mais avec des nombres flottants, après il faut gérer l’arrondi.

Exemple avec 5 lignes / 4 colonnes

Sommes de chaque ligne
30
22
24
19
27

Sommes de chaque colonne
36 26 31 29

Total 122

Valeurs calculées
8.852 6.393 7.623 7.131
6.492 4.689 5.590 5.230
7.082 5.115 6.098 5.705
5.607 4.049 4.828 4.516
7.967 5.754 6.861 6.418

EDIT: bon ça donne des valeurs plutôt "pondérées", si on veut que cela paraissent plus aléatoire on peut éventuellement appliquer un offset +/- à chaque valeur dont la somme est nulle en ligne / colonne.

fromvega

Si on prend la valeur arrondie, cela donne :

9 6 8 7
6 5 6 5
7 5 6 6
6 4 5 5
8 6 7 6

avec un total de 123, en regardant la somme des lignes et la somme des colonnes,
on identifie un élément à changer, ce qui conduit à la matrice :

9 6 8 7
6 5 6 5
7 5 6 6
6 4 4 5
8 6 7 6

Matrice que l’on eut ensuite bricoler avec la chaîne de Markov citée plus haut.

13 6 8 3
6 5 6 5
7 5 6 6
5 4 4 6
5 6 7 9

+0 -0

Un code Python à l’arrache, incomplet :
Je n’ai pas implémenté la modification de la matrice.
La différence entre le total cible et celui de la matrice
indique le nombre de cases à modifier.
Reste à identifier les cases en question, si il n’y en a qu’une, c’est facile.

### saisie contrainte lignes
print("contraintes sur les lignes :")
l = [0 for i in range(nbl)]
totl=0
for i in range(nbl):
    l[i]=int(input("ligne "+str(i)+":"))
    totl=totl+l[i]
print (l, "totl =", totl)

### saisie contraintes colonnes
c = [0 for i in range(nbc)]
totc=0
for i in range(nbc):
    c[i]=int(input("colonne "+str(i)+":"))
    totc=totc+c[i]
print (c, "totc =",totc)
if totc!=totl:
    print ("erreur, totc n'est pas égal à totl")

### calcul d'une matrice proche
for i in range(nbl):
    for j in range(nbc):                
        mat[i][j]=round(l[i]*c[j]/totl)
print(mat)
    
### calcul des ommes
lp = [0 for i in range(nbl)]
cp = [0 for i in range(nbc)]
tot=0
for i in range(nbl):
    for j in range(nbc):
        lp[i]=lp[i]+mat[i][j]
        cp[j]=cp[j]+mat[i][j]
        tot=tot+mat[i][j]
print (tot,lp, cp)
if tot!=totl:
    print("modier la matrice pour ",tot-totl," élément")
else:
    print ("la matrice est correcte")

Exemple

nombre de colonnes :2
nombre de lignes :3
contraintes sur les lignes :
ligne 0:21
ligne 1:17
ligne 2:35
[21, 17, 35] totl = 73
colonne 0:36
colonne 1:37
[36, 37] totc = 73
[[10, 11], [8, 9], [17, 18]]
73 [21, 17, 35] [35, 38]
la matrice est correcte
+0 -0

Matrice que l’on eut ensuite bricoler avec la chaîne de Markov citée plus haut.

Il doit y avoir un malentendu, la "chaîne de Markov" en question ne transforme pas n’importe quelle solution invalide en une solution valide — au contraire, elle est construite de manière à ne pas modifier les sommes sur les colonnes/lignes.


Comme c’est peut-être un peu abstrait, j’ai implémenté la méthode que j’ai décrite plus haut. Pour trouver une solution quelconque au problème, pas besoin de programmation linéaire : tout est dans la fonction initialize (vérifier par récurrence que ça marche toujours si la somme des contraintes sur les lignes et les colonnes est la même). Le 10000 est évidemment arbitraire — avis aux amateurs pour les meilleures conjectures.

Je précise que je ne suis pas Pythoniste de naissance. Ne copiez pas mon style ! ^^

import random

# Contraintes sur les lignes/colonnes
slin, scol = [30, 22, 24, 19, 27], [36, 26, 31, 29]
nb_lins, nb_cols = len(slin), len(scol)

# Génère deux entiers aléatoires distincts entre 0 et n-1 (uniformément)
def generate_distinct_pair(n):
    x = random.randrange(n)
    y = random.randrange(1, n)
    return x, (x + y) % n

# Effectue une étape de perturbation sur M
def perturbate(M):
    i1, i2 = generate_distinct_pair(nb_lins)
    j1, j2 = generate_distinct_pair(nb_cols)
    if M[i1][j1] != 0 and M[i2][j2] != 0:
        M[i1][j1] -= 1
        M[i2][j2] -= 1
        M[i1][j2] += 1
        M[i2][j1] += 1

# Retourne une matrice quelconque vérifiant les contraintes
def initialize():
    M = [[0] * nb_cols for _ in range(nb_lins)]
    for _ in range(nb_lins + nb_cols - 1):
        min_slin, min_scol = min(slin), min(scol)
        imin, jmin = slin.index(min_slin), scol.index(min_scol)
        if min_slin < min_scol:
            M[imin][jmin] = slin[imin]
            scol[jmin] -= slin[imin]
            slin[imin] = float('inf')
        else:
            M[imin][jmin] = scol[jmin]
            slin[imin] -= scol[jmin]
            scol[jmin] = float('inf')
    return M

M = initialize()
for _ in range(10000):
    perturbate(M)
print(M)

Matrice que l’on eut ensuite bricoler avec la chaîne de Markov citée plus haut.

Il doit y avoir un malentendu, la "chaîne de Markov" en question ne transforme pas n’importe quelle solution invalide en une solution valide — au contraire, elle est construite de manière à ne pas modifier les sommes sur les colonnes/lignes.

Non, pas de malentendu, je suis parti d’une matrice valide et j’ai généré une autre matrice valide.


Comme c’est peut-être un peu abstrait, j’ai implémenté la méthode que j’ai décrite plus haut. Pour trouver une solution quelconque au problème, pas besoin de programmation linéaire : tout est dans la fonction initialize (vérifier par récurrence que ça marche toujours si la somme des contraintes sur les lignes et les colonnes est la même). Le 10000 est évidemment arbitraire — avis aux amateurs pour les meilleures conjectures.

Mon code ne trouve pas toujours une solution valide :'( . Par contre, il part d’une idée assez intuitive ^^ .
Ton code trouve une solution valide :) , même si on a du mal à comprendre comment ça marche :magicien: .

Je précise que je ne suis pas Pythoniste de naissance. Ne copiez pas mon style ! ^^

Il ne faut pas non plus copier mon style. Pour ma part, je préfère le camelCase au snake-case que j’ai plus de mal à lire rapidement.

Source:Lucas-84

+1 -0

Non, pas de malentendu, je suis parti d’une matrice valide et j’ai généré une autre matrice valide.

Ah oui, au temps pour moi !

Ton code trouve une solution valide :) , même si on a du mal à comprendre comment ça marche :magicien: .

L’idée est en fait plutôt simple, quand on s’est débarrassé de ces contraintes ennuyantes "d’avoir l’air aléatoire" :


Au départ ta matrice est disons :

 ?  ?  ?  ?  30
 ?  ?  ?  ?  22
 ?  ?  ?  ?  24
 ?  ?  ?  ?  19
 ?  ?  ?  ?  27
36 26 31 29

On cherche la plus petite somme, que ce soit sur une ligne ou une colonne. Ici c’est 19, donc la quatrième ligne. On met un 19 n’importe où sur la ligne et 0 partout ailleurs, par exemple sur la première colonne. Sur la première colonne, la somme qu’on veut maintenant est 36–19 = 17.

 ?  ?  ?  ?  30
 ?  ?  ?  ?  22
 ?  ?  ?  ?  24
19  0  0  0   /
 ?  ?  ?  ?  27
17 26 31 29

Maintenant, on a le même le problème que précédemment, mais sur une matrice 4x4. Et on a conservé l’invariant que 17+26+31+29=30+22+24+27. La plus petite somme est 17 :

17  ?  ?  ?  13
 0  ?  ?  ?  22
 0  ?  ?  ?  24
19  0  0  0   /
 0  ?  ?  ?  27
 / 26 31 29

On a maintenant une matrice 4x3, et ainsi de suite. Les histoires de float('inf') c’est juste pour gérer les lignes/colonnes déjà remplies.


On notera au passage que cette stratégie donne une matrice de départ très "sparse" (de l’ordre de n+mn+m cases non nulles sur nmnm au total), alors que l’opération de perturbation ne progresse pas lorsqu’elle sélectionne certaines cases avec des zéros. Donc la partie initiale de la convergence n’est pas très optimisée.

L’idée est en fait plutôt simple, quand on s’est débarrassé de ces contraintes ennuyantes "d’avoir l’air aléatoire" :

Lucas-84

Merci, c’est très clair.
J’interprète la ligne du bas et la colonne de droite comme un "crédit" qu’il convient d’annuler.
La stratégie consiste à mettre à sec la ligne ou colonne la plus pauvre.
Cette stratégie permet de trouver une matrice valide en 9 itérations, elle est donc économe.
On peut envisager une stratégie (plus coûteuse en calcul) visant à appauvrir la case la plus riche, auquel cas on aurait 30 dans la case en haut à gauche.
30 est donc la plus grande valeur qui puisse apparaître dans cette matrice.

+0 -0
Connectez-vous pour pouvoir poster un message.
Connexion

Pas encore membre ?

Créez un compte en une minute pour profiter pleinement de toutes les fonctionnalités de Zeste de Savoir. Ici, tout est gratuit et sans publicité.
Créer un compte