Genera una matriz “aleatoria” de cierto rango sobre un conjunto fijo de elementos

Me gustaría generar matrices de tamaño m x n y rango r , con elementos provenientes de un conjunto finito específico, por ejemplo, {0,1} o {1,2,3,4,5} . Quiero que sean “aleatorios” en un sentido muy vago de esa palabra, es decir, quiero obtener una variedad de resultados posibles del algoritmo con una distribución vagamente similar a la distribución de todas las matrices sobre ese conjunto de elementos con el rango especificado.

De hecho, no me importa que tenga rango r , solo que está cerca de una matriz de rango r (medida por la norma de Frobenius).

Cuando el conjunto a la mano es el real, he estado haciendo lo siguiente, que es perfectamente adecuado para mis necesidades: generar matrices U de tamaño m x r y V de n x r , con elementos muestreados independientemente de, por ejemplo, Normal (0, 2). Entonces U V' es una matriz m x n de rango r (bueno, <= r , pero creo que es r con alta probabilidad).

Si solo hago eso y luego redondeo a binario / 1-5, el rango aumenta.

También es posible obtener una aproximación de rango inferior a una matriz haciendo una SVD y tomando los primeros r valores singulares. Sin embargo, esos valores no estarán en el conjunto deseado, y al redondearlos boostá nuevamente el rango.

Esta pregunta está relacionada, pero la respuesta aceptada no es “aleatoria”, y la otra respuesta sugiere SVD, que no funciona aquí como se señaló.

Una posibilidad que he pensado es hacer r vectores de fila o columna linealmente independientes del conjunto y luego obtener el rest de la matriz mediante combinaciones lineales de esos. Sin embargo, no estoy realmente claro, ni sobre cómo obtener vectores linealmente independientes lineales, o cómo combinarlos de forma casi aleatoria después de eso.

(No es que sea súper relevante, pero estoy haciendo esto en gran medida).


Actualización: He intentado el enfoque sugerido por EMS en los comentarios, con esta sencilla implementación:

 real = np.dot(np.random.normal(0, 1, (10, 3)), np.random.normal(0, 1, (3, 10))) bin = (real > .5).astype(int) rank = np.linalg.matrix_rank(bin) niter = 0 while rank > des_rank: cand_changes = np.zeros((21, 5)) for n in range(20): i, j = random.randrange(5), random.randrange(5) v = 1 - bin[i,j] x = bin.copy() x[i, j] = v x_rank = np.linalg.matrix_rank(x) cand_changes[n,:] = (i, j, v, x_rank, max((rank + 1e-4) - x_rank, 0)) cand_changes[-1,:] = (0, 0, bin[0,0], rank, 1e-4) cdf = np.cumsum(cand_changes[:,-1]) cdf /= cdf[-1] i, j, v, rank, score = cand_changes[np.searchsorted(cdf, random.random()), :] bin[i, j] = v niter += 1 if niter % 1000 == 0: print(niter, rank) 

Funciona rápidamente para matrices pequeñas, pero se desmorona para, por ejemplo, 10×10; parece que se queda atascado en el rango 6 o 7, al menos en cientos de miles de iteraciones.

Parece que esto podría funcionar mejor con una función objective mejor (es decir, menos plana), pero no sé qué sería eso.


También he intentado un método de rechazo simple para construir la matriz:

 def fill_matrix(m, n, r, vals): assert m >= r and n >= r trans = False if m > n: # more columns than rows I think is better m, n = n, m trans = True get_vec = lambda: np.array([random.choice(vals) for i in range(n)]) vecs = [] n_rejects = 0 # fill in r linearly independent rows while len(vecs)  len(vecs): vecs.append(v) else: n_rejects += 1 print("have {} independent ({} rejects)".format(r, n_rejects)) # fill in the rest of the dependent rows while len(vecs)  len(vecs): n_rejects += 1 if n_rejects % 1000 == 0: print(n_rejects) else: vecs.append(v) print("done ({} total rejects)".format(n_rejects)) m = np.vstack(vecs) return mT if trans else m 

Esto funciona bien, por ejemplo, para matrices binarias 10×10 con cualquier rango, pero no para matrices 0-4 o binarios mucho más grandes con rango inferior. (Por ejemplo, obtener una matriz binaria de 20×20 de rango 15 me llevó 42,000 rechazos; con 20×20 de rango 10, tomó 1.2 millones).

Esto se debe claramente a que el espacio que abarcan las primeras r filas es una porción muy pequeña del espacio del que tomo muestras, por ejemplo, {0,1}^10 , en estos casos.

Queremos la intersección del tramo de las primeras r filas con el conjunto de valores válidos. Así que podríamos intentar muestrear el intervalo y buscar valores válidos, pero dado que el intervalo implica coeficientes de valores reales, nunca nos encontrarán vectores válidos (incluso si lo normalizamos de modo que, por ejemplo, el primer componente esté en el conjunto válido).

Tal vez esto puede ser formulado como un problema de progtwigción de enteros, o algo así?

A mi amigo, Daniel Johnson, que comentó anteriormente, se le ocurrió una idea, pero veo que nunca la publicó. No es muy sencillo, pero es posible que puedas adaptarlo.

Si A es m-by-r y B es r-by-n y ambos tienen rango r, entonces AB tiene rango r. Ahora, solo tenemos que elegir A y B para que AB tenga valores solo en el conjunto dado. El caso más simple es S = {0,1,2,...,j} .

Una opción sería hacer A binario con sums de fila / col apropiadas que garanticen el rango correcto y B con sums de columna que no j más de j (de modo que cada término en el producto esté en S ) y sums de fila elegidas para causar el rango r (o al menos alentarlo ya que puede usarse el rechazo).

Solo creo que podemos llegar a dos esquemas de muestreo independientes en A y B que son menos complicados y más rápidos que intentar atacar toda la matriz a la vez. Desafortunadamente, todo mi código de muestreo matricial está en la otra computadora. Sé que se generaliza fácilmente para permitir entradas en un conjunto mayor que {0,1} (es decir, S ), pero no puedo recordar cómo se escalaron los cálculos con m*n .

¿Qué tal esto?

 rank = 30 n1 = 100; n2 = 100 from sklearn.decomposition import NMF model = NMF(n_components=rank, init='random', random_state=0) U = model.fit_transform(np.random.randint(1, 5, size=(n1, n2))) V = model.components_ M = np.around(U) @ np.around(V) 

No estoy seguro de cuán útil será esta solución, pero puede construir una matriz que le permita buscar la solución en otra matriz con solo 0 y 1 como entradas. Si busca aleatoriamente en la matriz binaria, es equivalente a modificar aleatoriamente los elementos de la matriz final, pero es posible encontrar algunas reglas para hacerlo mejor que una búsqueda aleatoria.

Si desea generar una matriz m -by- sobre el conjunto de elementos E con los elementos e i , 0<=i , comience con la matriz m -by- k*m , A :

Matriz generadora

Claramente, esta matriz tiene rango m . Ahora, puedes construir otra matriz, B, que tenga 1s en ciertas ubicaciones para elegir los elementos del conjunto E. La estructura de esta matriz es:

Matriz selectora

Cada B i es una matriz k by- n . Entonces, el tamaño de A B es m -by- n y el rango ( A B ) es min (m, rango ( B )) . Si queremos que la matriz de salida tenga solo elementos de nuestro conjunto, E , entonces cada columna de B i debe tener exactamente un elemento establecido en 1 , y el rest establecido en 0 .

Si desea buscar un determinado rango en B al azar, debe comenzar con un B válido con rango máximo y rotar una columna aleatoria j de un B i aleatorio por una cantidad aleatoria. Esto es equivalente a cambiar la columna i fila j de A * B a un elemento aleatorio de nuestro conjunto, por lo que no es un método muy útil.

Sin embargo, puedes hacer ciertos trucos con las matrices. Por ejemplo, si k es 2 y no hay superposiciones en las primeras filas de B 0 y B 1 , puede generar una fila dependiente linealmente agregando las primeras filas de estas dos submatrices. La segunda fila también dependerá linealmente de las filas de estas dos matrices. No estoy seguro de si esto se generalizará fácilmente a k más grande que 2, pero estoy seguro de que habrá otros trucos que puede emplear.

Por ejemplo, un método simple para generar a lo sumo el rango k (cuando m es k+1 ) es obtener una B 0 válida al azar, seguir girando todas las filas de esta matriz hasta obtener B 1 a B m-2 , configurar la primera fila de B m-1 a todos 1, y las filas restantes a todos 0. El rango no puede ser menor que k (suponiendo que n > k ), porque las columnas B_0 tienen exactamente 1 elemento distinto de cero. Las filas restantes de las matrices son todas combinaciones lineales (de hecho, copias exactas para casi todas las submatrices) de estas filas. La primera fila de la última submatriz es la sum de todas las filas de la primera submatriz, y las filas restantes de la misma son todas ceros. Para valores mayores de m , puede usar permutaciones de filas de B 0 en lugar de una simple rotación.

Una vez que genere una matriz que satisfaga la restricción de rango, puede evitar mezclar aleatoriamente las filas y columnas de la misma para generar otras.