Сгенерировать «случайную» матрицу определенного ранга по фиксированному набору элементов

Я хотел бы сгенерировать матрицы размера mx nи ранга rс элементами, поступающими из заданный конечный набор, например {0,1}или {1,2,3,4,5}. Я хочу, чтобы они были «случайными» в каком-то очень широком смысле этого слова, т. Е. Я хочу получить множество возможных выходов из алгоритма с распределением, отдаленно похожим на распределение всех матриц по этому набору элементов с указанным рангом.

На самом деле меня не волнует, что она имеет ранг r, просто она близка к матрице ранга r(, измеряемой нормой Фробениуса ).

Когда под рукой имеется множество вещественных чисел, я делаю следующее, что вполне соответствует моим потребностям :генерирую матрицы Uразмера mx rи Vиз nx r, с элементами, взятыми независимо, например, из Нормальный (0, 2 ). Тогда U V'— это матрица mx nрангаr(ну, <= r, но я думаю, что rс высокой вероятностью ).

Если я просто сделаю это, а затем округлю до двоичного / 1 -5, тем не менее, ранг увеличится.

Также можно получить приближение более низкого -ранга к матрице, выполнив SVD и взяв первые rсингулярные значения. Эти значения, однако, не будут лежать в желаемом наборе, и их округление снова повысит ранг.

Этот вопрос связан, но принятый ответ не является «случайным», а другой ответ предполагает SVD, который, как уже отмечалось, здесь не работает.

Одна возможность, о которой я подумал, состоит в том, чтобы сделать rлинейно независимыми векторы-строки или столбцы из набора, а затем получить остальную часть матрицы линейными комбинациями этих. Однако я не очень понимаю, как получить «случайные» линейно независимые векторы или как после этого комбинировать их квазислучайным образом.

(Не то, чтобы это супер -актуально, но я делаю это в numpy.)

Обновление:Я попробовал подход, предложенный EMS в комментариях, с этой простой реализацией:

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)

Он работает быстро для небольших матриц, но разваливается, например, для. 10x10 --кажется, что он застревает на 6 или 7 ранге, по крайней мере, на сотни тысяч итераций.

Кажется, что это могло бы работать лучше с лучшей (т.е. менее -плоской )целевой функцией, но я не знаю, что это будет.

Я также попробовал простой метод отбраковки для построения матрицы.:

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) < r:
        v = get_vec()
        if np.linalg.matrix_rank(np.vstack(vecs + [v])) > 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) < m:
        v = get_vec()
        if np.linalg.matrix_rank(np.vstack(vecs + [v])) > 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 m.T if trans else m

Это хорошо работает, например, для. Двоичные матрицы 10x10 с любым рангом, но не для 0 -4 матриц или гораздо больших двоичных матриц с более низким рангом. (Например, для получения бинарной матрицы 20x20 ранга 15 мне потребовалось 42 000 отказов; с 20х20 10 ранга ушло 1,2 млн.)

Очевидно, это связано с тем, что пространство, занимаемое первыми rстроками, является слишком малой частью пространства, из которого я делаю выборку, например. {0,1}^10, в этих случаях.

Нам нужно пересечение диапазона первых rстрок с набором допустимых значений. Таким образом, мы могли бы попробовать выполнить выборку из диапазона и найти допустимые значения,но поскольку диапазон включает в себя реальные -значения коэффициентов, которые никогда не найдут нам действительные векторы (, даже если мы нормализуем так, что, например. первый компонент находится в допустимом наборе ).

Может быть, это можно сформулировать как задачу целочисленного программирования или что-то в этом роде?