Delete rows or columns of matrix containing invalid elements, such that a maximum number of valid elements is kept
Here is a solution based on back-tracking. You probably need to double check but it seemed ok on my side. Should take 1h on a 30x900 matrix, but you can use multi-threading in the bottom for loop to speed it up.
import numpy as np
from tqdm import tqdm


def hash_list(l1, l2):
    return str(sorted(l1)) + str(sorted(l2))  # could be faster


def max_matrix(a):
    # the process:
    # we construct list of rows and columns to select, which keep only valid values (a[rows, columns] = 1 everywhere)
    # at each step, for each other row/column, if it can be added in a valid way, we carry on deeper
    # when it's not possible to extend anymore
    # we backtrack to where we were and try to add the next possible row/column.
    # we store in res the visited valid rows/columns, which prevents most redundant computations
    res = {}
    M, N = a.shape

    def make_grow(a, sel_columns, sel_rows, other_columns, other_rows, d='right'):
        assert len(sel_columns) + len(other_columns) == M
        assert len(sel_rows) + len(other_rows) == N
        h = hash_list(sel_columns, sel_rows) + d
        if h in res: # we went through this already
            return

        # attempting to add a column
        if d == 'right':
            for i in other_columns:
                if 0 not in a[i, sel_rows]:
                    # we can append it: let's try to go left or right from there
                    other_columns.remove(i)
                    sel_columns.append(i)
                    make_grow(a, sel_columns, sel_rows, other_columns, other_rows, d='left')
                    make_grow(a, sel_columns, sel_rows, other_columns, other_rows, d='right')
                    
        # attempting to add a row
        else:
            for j in other_rows:
                if 0 not in a[sel_columns, j]:
                    # we can append it: let's try to go left or right from there
                    other_rows.remove(j)
                    sel_rows.append(j)
                    make_grow(a, sel_columns, sel_rows, other_columns, other_rows, d='right')
                    make_grow(a, sel_columns, sel_rows, other_columns, other_rows, d='left')

        # this position has been visited.
        res[h] = (sel_columns, sel_rows)

    # now launching from all possible first positions.
    for i in tqdm(range(M)):
        for j in range(N):
            if a[i, j] == 1:
                make_grow(a, [i], [j], [elt for elt in range(M) if elt != i], [elt for elt in range(N) if elt != j], d='right')
                make_grow(a, [i], [j], [elt for elt in range(M) if elt != i], [elt for elt in range(N) if elt != j], d='left')
            
    return res


if __name__ == '__main__':
    a = np.random.uniform(size=(30, 900))
    a[a &gt; 0.3] = 1
    a[a &lt;= 0.3] = 0
    a = a.astype(np.uint8)
    print(a)
    res = max_matrix(a)

    # just looking through the results to see the best configuration. (all valid configurations are stored btw)
    res_max = 0
    l1_max, l2_max = [], []

    for _, (l1, l2) in res.items():
        if len(l1) * len(l2) &gt; res_max:
            res_max = len(l1) * len(l2)
            l1_max = l1
            l2_max = l2
    print(res_max, l1_max, l2_max)