In general, I don't think you can expect a set to be well approximated by such a small number of rectangles.

Let $S$ be a random set formed by including every square with probability $1/2$. Then with high probability $S$ has $r \geq 0.5-\epsilon$ for any $\epsilon$.

Now consider any (fixed) arbitrary set $T$. The error of $T$ from $S$ can be thought of as the sum of $N^2$ Bernoulli trials, each with probability $1/2$. It follows from the Chernoff bound that for any fixed $\epsilon<1/2$ the probability of having error at most $\epsilon N^2$ is at most $c^{n^2}$ for some constant $c<1$ depending only on $\epsilon$.

On the other hand, there are only $4^N$ rectangles (choose whether or not to include each row and column), so at most $4^{Nk}$ unions of $k$ rectangles. Taking the union bound over all such unions, we see that with high probability $S$ is not approximated by any union of $o(N)$ rectangles.

In general, I have a feeling (though I'm not familiar enough with this area to say for certain) that a better way to explain this all is in terms of information theory/entropy -- Specifying that a set has density $r$ still leaves you with approximately $N^2 H(r)$ (where $H$ is the entropy function) bits of entropy. On the other hand, the union of $k$ rectangles has less than $2kN$ such bits. You can't compress the former into the latter if $k$ is much less than $2kH(r)$ without incurring a fair amount of error.

In general, I don't think you can expect a set to be well approximated by such a small number of rectangles.

Let $S$ be a random set formed by including every square with probability $1/2$. Then with high probability $S$ has $r \geq 0.5-\epsilon$ for any $\epsilon$.

Now consider any (fixed) arbitrary set $T$. The error of $T$ from $S$ can be thought of as the sum of $N^2$ Bernoulli trials, each with probability $1/2$. It follows from the Chernoff bound that for any fixed $\epsilon<1/2$ the probability of having error at most $\epsilon N^2$ is at most $c^{n^2}$ for some constant $c<1$ depending only on $\epsilon$.

On the other hand, there are only $4^N$ rectangles (choose whether or not to include each row and column), so at most $4^{Nk}$ unions of $k$ rectangles. Taking the union bound over all such unions, we see that with high probability $S$ is not approximated by any union of $o(N)$ rectangles.

In general, I have a feeling (though I'm not familiar enough with this area to say for certain) that a better way to explain this all is in terms of information theory/entropy -- Specifying that a set has density $r$ still leaves you with approximately $N^2 H(r)$ (where $H$ is the entropy function) bits of entropy. On the other hand, the union of $k$ rectangles has less than $2kN$ such bits. You can't compress the former into the latter if $k$ is much less than $NH(r)$ without incurring a fair amount of error.