在Perl中随机化矩阵,同时保持行和列总和不变

11
我有一个矩阵,我想要随机它几千次,同时保持行和列的总数不变:
     1 2 3 
   A 0 0 1 
   B 1 1 0 
   C 1 0 0      

一个有效的随机矩阵示例如下:

     1 2 3
   A 1 0 0
   B 1 1 0
   C 0 0 1

我的实际矩阵规模要大得多(大约有600x600个项目),因此我真的需要一种计算效率高的方法。

我的最初(低效)的方法是使用Perl Cookbook中的shuffle函数来打乱数组。

我在下面粘贴了我的当前代码。如果在while循环中找不到解决方案,我已经加入了额外的代码以从新的随机数字列表开始。这个算法对于小矩阵来说运行得很好,但是一旦我开始扩展它,它就需要很长时间才能找到符合要求的随机矩阵。

有没有更有效的方法来完成我正在寻找的东西? 非常感谢!

#!/usr/bin/perl -w
use strict;

my %matrix = ( 'A' => {'3'  => 1 },
           'B' => {'1'  => 1,
               '2'  => 1 },
           'C' => {'1'  => 1 }
    );

my @letters = ();
my @numbers = ();

foreach my $letter (keys %matrix){
    foreach my $number (keys %{$matrix{$letter}}){
    push (@letters, $letter);
    push (@numbers, $number);
    }
}

my %random_matrix = ();

&shuffle(\@numbers);
foreach my $letter (@letters){
    while (exists($random_matrix{$letter}{$numbers[0]})){
    &shuffle (\@numbers);
    }
    my $chosen_number = shift (@numbers);
    $random_matrix{$letter}{$chosen_number} = 1;
}

sub shuffle {
    my $array = shift;
    my $i = scalar(@$array);
    my $j;
    foreach my $item (@$array )
    {
        --$i;
        $j = int rand ($i+1);
        next if $i == $j;
        @$array [$i,$j] = @$array[$j,$i];
    }
    return @$array;
}

2
是的,矩阵始终是二进制的。此外,0的实例数量远远超过1的实例数量! - Lucas
3
请勿发布未经授权的书籍副本链接!已被标记为垃圾信息。 - Sinan Ünür
2
@Lucas:在这个过程中,你覆盖了我在Perl FAQ列表中关于shuffle的链接,只是用一个指向某个网页托管公司的404链接替换了它。我会说你在这里继续发送垃圾信息。 - Sinan Ünür
5
@Sinan + Lucas:即使没有书籍链接,这个问题仍然有意义,对吗?我相信大多数人都能理解洗牌是怎么回事,如果不行,还可以用谷歌搜索。这样就不用争论哪个链接最好了 :) - Mark Byers
2
@Sinan,我非常抱歉,这是我在这里的第一个问题,所以a)我没有意识到您编辑了我的链接(甚至不知道会发生这种情况),b)我试图删除我的链接到非法网站,并用正确的链接替换它。感谢您恢复它的努力,但是称我为“垃圾邮件”可能有点太快了! - Lucas
显示剩余7条评论
6个回答

10
你当前算法的问题在于,你试图通过洗牌来摆脱死胡同——具体来说,当你的@letters@numbers数组(在@numbers的初始洗牌后)产生多次相同的单元格时。这种方法在矩阵很小的情况下有效,因为不需要尝试太多次才能找到可行的重新洗牌。然而,在列表很大的情况下,这种方法是致命的。即使你可以更高效地寻找替代方案——例如,尝试排列而不是随机洗牌——这种方法也可能注定失败。
与其洗牌整个列表,不如通过对现有矩阵进行小的修改来解决问题。
例如,让我们从你的示例矩阵(称之为M1)开始。随机选择一个单元格进行更改(比如A1)。此时,矩阵处于非法状态。我们的目标是以最少的编辑次数来修复它——具体来说是3次额外的编辑。你可以通过“走动”矩阵来实现这3次额外的编辑,每次修复一行或一列都会产生另一个需要解决的问题,直到你走完整个矩形。
例如,在将A1从0更改为1后,有3种方法可以进行下一次修复:A3、B1和C1。让我们决定第一个编辑应该修复行。所以我们选择A3。在第二个编辑中,我们将修复列,所以我们有选择:B3或C3(假设是C3)。最后的修复只提供了一个选择(C1),因为我们需要返回到原始编辑的列。最终结果是一个新的有效矩阵。
    Orig         Change A1     Change A3     Change C3     Change C1
    M1                                                     M2

    1 2 3        1 2 3         1 2 3         1 2 3         1 2 3
    -----        -----         -----         -----         -----
A | 0 0 1        1 0 1         1 0 0         1 0 0         1 0 0
B | 1 1 0        1 1 0         1 1 0         1 1 0         1 1 0
C | 1 0 0        1 0 0         1 0 0         1 0 1         0 0 1

如果编辑路径导致死胡同,你需要回溯。如果所有修复路径都失败,最初的编辑可能会被拒绝。

这种方法可以快速生成新的有效矩阵。但不一定会产生随机结果:M1和M2仍然高度相关,当矩阵大小增长时,这一点将变得更加明显。

如何增加随机性?你提到大多数单元格(99%或更多)都是零。一个想法是按照以下步骤进行:对于矩阵中的每个1,将其值设置为0,然后使用上面概述的4次编辑方法来修复矩阵。实际上,你将把所有1移到新的随机位置。

这里有一个例子。可能还有进一步的速度优化,但这种方法在我的Windows电脑上,在30秒左右内以0.5%的密度产生了10个新的600x600矩阵。不知道是否足够快。

use strict;
use warnings;

# Args: N rows, N columns, density, N iterations.
main(@ARGV);

sub main {
    my $n_iter = pop;
    my $matrix = init_matrix(@_);
    print_matrix($matrix);
    for my $n (1 .. $n_iter){
        warn $n, "\n"; # Show progress.
        edit_matrix($matrix);
        print_matrix($matrix);
    }
}

sub init_matrix {
    # Generate initial matrix, given N of rows, N of cols, and density.
    my ($rows, $cols, $density) = @_;
    my @matrix;
    for my $r (1 .. $rows){
        push @matrix, [ map { rand() < $density ? 1 : 0  } 1 .. $cols ];
    }
    return \@matrix;
}

sub print_matrix {
    # Dump out a matrix for checking.
    my $matrix = shift;
    print "\n";
    for my $row (@$matrix){
        my @vals = map { $_ ? 1 : ''} @$row;
        print join("\t", @vals), "\n";
    }
}

sub edit_matrix {
    # Takes a matrix and moves all of the non-empty cells somewhere else.
    my $matrix = shift;
    my $move_these = cells_to_move($matrix);
    for my $cell (@$move_these){
        my ($i, $j) = @$cell;
        # Move the cell, provided that the cell hasn't been moved
        # already and the subsequent edits don't lead to a dead end.
        $matrix->[$i][$j] = 0
            if $matrix->[$i][$j]
            and other_edits($matrix, $cell, 0, $j);
    }
}

sub cells_to_move {
    # Returns a list of non-empty cells.
    my $matrix = shift;
    my $i = -1;
    my @cells = ();
    for my $row (@$matrix){
        $i ++;
        for my $j (0 .. @$row - 1){
            push @cells, [$i, $j] if $matrix->[$i][$j];
        }
    }
    return \@cells;
}

sub other_edits {
    my ($matrix, $cell, $step, $last_j) = @_;

    # We have succeeded if we've already made 3 edits.
    $step ++;
    return 1 if $step > 3;

    # Determine the roster of next edits to fix the row or
    # column total upset by our prior edit.
    my ($i, $j) = @$cell;
    my @fixes;
    if ($step == 1){
        @fixes = 
            map  { [$i, $_] }
            grep { $_ != $j and not $matrix->[$i][$_] }
            0 .. @{$matrix->[0]} - 1
        ;
        shuffle(\@fixes);
    }
    elsif ($step == 2) {
        @fixes = 
            map  { [$_, $j] }
            grep { $_ != $i and $matrix->[$_][$j] }
            0 .. @$matrix - 1
        ;
        shuffle(\@fixes);
    }
    else {
        # On the last edit, the column of the fix must be
        # the same as the column of the initial edit.
        @fixes = ([$i, $last_j]) unless $matrix->[$i][$last_j];
    }

    for my $f (@fixes){
        # If all subsequent fixes succeed, we are golden: make
        # the current fix and return true.
        if ( other_edits($matrix, [@$f], $step, $last_j) ){
            $matrix->[$f->[0]][$f->[1]] = $step == 2 ? 0 : 1;
            return 1;
        }
    }

    # Failure if we get here.
    return;
}

sub shuffle {
    my $array = shift;
    my $i = scalar(@$array);
    my $j;
    for (@$array ){
        $i --;
        $j = int rand($i + 1);
        @$array[$i, $j] = @$array[$j, $i] unless $i == $j;
    }
}

1
算法绝对够快!但我有点担心它的非随机性,考虑到限制条件,我需要尽可能随机的矩阵。我会尝试按照你的算法生成矩阵,并检查它们与原始矩阵的相关性。问题是,当然我无法将它们与真正随机的矩阵进行比较 ;)。 - Lucas
1
@Lucas 算法在代码中的实现非常随机。原始矩阵中的每个1都被关闭,另外3个单元格被切换以保持行/列总和不变。这3个修复步骤尽可能随机,但仍需满足问题的约束条件。如果有一种非随机的方面,那就是因为原始矩阵中的1在新矩阵中变成0的概率高于随机。纠正这种偏差的一种方法可能是跳过一些比例的原始矩阵中的1 - 换句话说,将它们保留为原值。这个比例是多少?也许可以使用密度来计算。 - FMc
@Lucas:如果您选择此方法,改善随机性的另一种方法是多次重复随机化过程(随机次数)。然后我认为与原矩阵的相关性变得可以忽略不计。这会使它稍微慢一些,但实现比我的建议更简单。 - Mark Byers
是的,昨晚我在各种随机矩阵上运行了这个算法几千次,并看到每次迭代相关性都稳步下降。 很好的解决方案FM! 您可以解释一下如何在init矩阵中使用问号映射函数吗?您正在生成一个数组的数组作为矩阵,但我对此符号不熟悉。非常感谢您为我的问题提供了一个很好的解决方案:)。 - Lucas
@Lucas 很高兴能帮上忙。这是一个有趣的问题,与我最近正在进行的一个副项目有些相关。关于 map { rand() < $density ? 1 : 0 } 1 .. $cols,我们正在生成一个由 0 和 1 组成的列表。如果 rand() 小于 $density,则该单元格将为 1;否则为 0。搜索 perldoc perlop 以获取有关 ? : 语法的条件运算符的更多详细信息。它基本上是一个迷你 IF-THEN 结构,可以方便地在较大的表达式中使用。 - FMc

6
步骤1:首先,我会将矩阵初始化为零,并计算所需的行和列总数。
步骤2:现在选择一个随机行,其权重由必须在该行中包含的1的数量决定(因此,具有300个计数的行比具有5个计数的行更有可能被选择)。
步骤3:对于此行,选择一个随机列,其权重由该列中包含的1的数量决定(除了忽略可能已经包含1的单元格-稍后再说)。
步骤4:将1放置在此单元格中,并减少适当行和列的行和列计数。
步骤5:返回步骤2,直到没有行具有非零计数。
问题在于,此算法可能无法终止,因为您可能有需要放置1的行和需要1的列,但您已经在该单元格中放置了1,因此您会“卡住”。我不确定这种情况发生的可能性有多大,但我不会感到惊讶,如果它发生得非常频繁-足以使算法无法使用。如果这是一个问题,我可以想出两种解决方法:
a)递归构造上述算法,并允许在失败时回溯。
b)如果没有其他选择,则允许单元格包含大于1的值并继续。然后,在最后,您具有正确的行和列计数,但某些单元格可能包含大于1的数字。您可以通过找到类似于此的分组来解决此问题:
2 . . . . 0
. . . . . .
. . . . . .
0 . . . . 1

并将其更改为:

1 . . . . 1
. . . . . .
. . . . . .
1 . . . . 0

如果您有很多零,那么找到这样的分组应该很容易。我认为b)可能更快。

我不确定这是否是最好的方法,但它可能比打乱数组要快。我将跟踪这个问题,看看其他人想出了什么。


谢谢您的建议!这似乎是一个可行的方法,我一定会尝试实现。我将尝试查看算法需要填充已经填充的单元格的次数。 由于我的表中1的频率很低,我不指望这种情况经常发生,所以再次运行算法也可能是可行的(如果它不经常发生)! - Lucas
如果使用算法b),7并不是问题。您只需要重复足够多次的重新组织,每次选择不同的块即可。 - Mark Byers
在选择随机列时,请记得不要包括已经包含1的列。 - Mark Byers
(除非您别无选择) - Mark Byers
我没有对随机选择进行加权,所以现在我要实现它。今晚我会在我的网站上发布一个样本矩阵,当它准备好时,我会在这里发布链接。 - Lucas
显示剩余2条评论

1

我不是数学家,但我认为如果您需要保持相同的列和行总数,则矩阵的随机版本将具有相同数量的1和0。

如果我错了,请纠正我,但这意味着制作矩阵的后续版本只需要您重新排列行和列即可。

随机洗牌列不会改变行和列的总数,随机洗牌行也不会。因此,我会先洗牌行,然后再洗牌列。

这应该非常快。


1
我的理解是行(和列)的总数必须保持恒定且顺序不变。如果是这样,那么交换两行是可以的,但是只有当它们具有相同的行总数时才可以,如果它们具有不同的行总数则不可以。在问题中给出的 3x3 示例中,交换 A 和 C 行是有效的,因为它们都具有行总数 1,但是不能交换 A 和 B 行,因为行 B 的总数为 2。 - Mark Byers
随机洗牌列不会影响行总数,但会影响列总数(反之亦然)。如果将初始矩阵中的所有列向右移动一位,则第1列的总和不再等于2。 编辑:@Mark 您是正确的,而且您的例子表述得比我好 :) - Lucas
1
此外,有些可能需要生成的解决方案并不能仅通过行和列交换来实现。例如,如果您有一个4x4矩阵,其中行和列的总数均为2,并且初始配置为1100,1100,0011,0011,我不确定您如何仅通过行和列交换将其更改为1100,1010,0101,0011 - Mark Byers
是的,我认为问题描述有点混淆。我最初的印象是需要保持每列和每行的总和不变,而不是总和的位置。经过另一次查看,随机洗牌列和行是行不通的。 - Tim Rupe

0
不确定它是否有帮助,但你可以尝试从一个角落开始,每列和每行都要跟踪总数和实际总数。不要试图找到完美的矩阵,而是将总数视为金额并拆分它。对于每个元素,找到行总数-实际行总数和列总数-实际列总数中较小的数字。现在你有了随机数的上限。 清楚吗?很抱歉我不懂Perl,所以不能展示任何代码。

你的回答不是很清楚,但我现在明白了。我认为你错过了其中一个要求:每个单元格只能包含0或1。虽然问题中没有明确说明,但在问题的评论中提到了这一点。 - Mark Byers

0

像 @Gabriel 一样,我不是 Perl 程序员,所以你的代码可能已经实现了这个功能...

你只发布了一个示例。不清楚你想要一个随机矩阵,该矩阵每行和每列都有相同数量的 1,还是一个具有相同行和列但被打乱的矩阵。如果后者足够好,你可以创建一个行(或列,无所谓)索引数组,并随机排列它。然后,你可以按照随机化的索引指定的顺序读取原始数组。无需修改原始数组或创建副本。

当然,这可能无法满足你未明确说明的要求方面。


谢谢Mark!但是正如你可能已经在其他评论中读到的那样,我正在寻找解决你所描述的前一个问题的解决方案 :)。 - Lucas

0
感谢FMc的Perl代码。基于这个解决方案,我将其重写为Python(供自己使用并在此分享以获得更清晰的表述),如下所示:
matrix = numpy.array( 
    [[0, 0, 1], 
     [1, 1, 0], 
     [1, 0, 0]]
)

def shuffle(array):
    i = len(array)
    j = 0
    for _ in (array):
        i -= 1;
        j = random.randrange(0, i+1) #int rand($i + 1);
        #print('arrary:', array)
        #print(f'len(array)={len(array)}, (i, j)=({i}, {j})')
        if i != j: 
            tmp = array[i]
            array[i] = array[j]
            array[j] = tmp
    return array

def other_edits(matrix, cell, step, last_j):
    # We have succeeded if we've already made 3 edits.
    step += 1
    if step > 3: 
        return True

    # Determine the roster of next edits to fix the row or
    # column total upset by our prior edit.
    (i, j) = cell
    fixes = []
    if (step == 1):
        fixes = [[i, x] for x in range(len(matrix[0])) if x != j and not matrix[i][x] ]
        fixes = shuffle(fixes)
    elif (step == 2):
        fixes = [[x, j] for x in range(len(matrix)) if x != i and matrix[x][j]]
        fixes = shuffle(fixes)
    else:
        # On the last edit, the column of the fix must be
        # the same as the column of the initial edit.
        if not matrix[i][last_j]: fixes = [[i, last_j]]

    for f in (fixes):
        # If all subsequent fixes succeed, we are golden: make
        # the current fix and return true.
        if ( other_edits(matrix, f, step, last_j) ):
            matrix[f[0]][f[1]] = 0 if step == 2 else 1
            return True

    # Failure if we get here.
    return False # return False

def cells_to_move(matrix):
    # Returns a list of non-empty cells.
    i = -1
    cells = []
    for row in matrix:
        i += 1;
        for j in range(len(row)):
            if matrix[i][j]: cells.append([i, j])
    return cells

def edit_matrix(matrix):
    # Takes a matrix and moves all of the non-empty cells somewhere else.
    move_these = cells_to_move(matrix)
    for cell in move_these:
        (i, j) = cell
        # Move the cell, provided that the cell hasn't been moved
        # already and the subsequent edits don't lead to a dead end.
        if matrix[i][j] and other_edits(matrix, cell, 0, j):
            matrix[i][j] = 0
    return matrix

def Shuffle_Matrix(matrix, N, M, n_iter):
    for n in range(n_iter):
        print(f'iteration: {n+1}') # Show progress.
        matrix = edit_matrix(matrix)
        #print('matrix:\n', matrix)
    return matrix

print(matrix.shape[0], matrix.shape[1]) 

# Args: N rows, N columns, N iterations.
matrix2 = Shuffle_Matrix(matrix, matrix.shape[0], matrix.shape[1], 1) 

print("The resulting matrix:\n", matrix2)

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接