Generative model that satisfies certain algebraic constraints

Question

Disclaimer: I need guidance and help with where to start looking for solution of the problem I have described below. My background is in optimization and I am new to statistical methods, so there is a good chance that I am asking the wrong question or/and used wrong terminologies (please correct if that is the case).

Below I setup my problem:

Given:

• the set of $n\times n$ matrices.

• two functoins, $f: \mathbb{R}^{n^2} \rightarrow\mathbb{R}$ and ${\bf{g}}: \mathbb{R}^{n^2}\rightarrow Symm.(n\times n)$

• Two constraints as follows:

$$1 - \epsilon_1<f({\bf{M}}) < 1 + \epsilon_1$$

$$-\epsilon_2< \lvert\lvert{{\bf{g}}}({\bf{M}})\rvert\rvert_{\infty} < \epsilon_2$$

for ${\bf{M}} \in \mathbb{R}^{n^2}$ and $\epsilon_1$ and $\epsilon_2$ both fixed small positive numbers.

Here are my questions:

1. Can I find a model, or a distribution, which when I sample from it, it produces $n\times n$ matrices that satisfy the above two constraints (most of the time)? The sampled data needs to be close the real distribution in order to be representative.
2. Is question (1) a well formulated question?
3. If the answer to (2) is yes, what method(s) should I look into in order to work towards a solution?

For those who are interested in more concrete realizations of functions $f$ and $\bf{g}$, $f({\bf{M}})=\mathrm{det}({\bf{M}})$ and ${{\bf{g}}}({\bf{M}})=\frac{1}{2}({\bf{M}}^T{\bf{M}} - {\bf{I}})$.

I appreciate any hint or help with this problem. Tags mentioned below are speculative.

Answer 1

Your question defines a set of matrices you want to sample from. However as pointed out in the comments, the notion of sampling depends on some distribution or probability measure over the set.

As a simple example, consider the $\delta$ distribution around some matrix $M_0$ that satisfies your constraints - this is a perfectly valid distribution and it answers your requirements because 'sampling' from it will always give you the matrix $M_0$, but this is probably not what you want.

Now having that in mind we can try to come up with of with some distribution that might be reasonable for your purpose (with the caveat that without knowing the purpose of this exercise it is hard to tell what is 'reasonable'): note that your specific set of constraints defines matrices that are close to unitary, since if $M$ is unitary ($M^TM=I$) then both $f(M)=1$ and $g(M)=0$. So, let

$$M = U(I + \epsilon X)$$

Where $U$ is unitary, $X$ is a random matrix (with some distribution) and $\epsilon$ is a small real number, therefore

$$ \det(M)=1 + \epsilon \text{tr}(X) + O(\epsilon^2) $$

$$M^TM - I = \epsilon(X+X^T) + O(\epsilon^2) $$

For a given $X$ this defines a possible range for $\epsilon$, since you require that

$$ -\epsilon_1 < \epsilon \text{tr}(X) < \epsilon_2 $$ $$ -\epsilon_1 < |\epsilon| \cdot \|\frac{X+X^T}{2}\|_{\infty} < \epsilon_2 $$

(For the second constraints only the upper bound is actually relevant since the norm can only be positive). So, generate a random matrix $X$ form some distribution, then generate a random $\epsilon$ from the range defined above (or alternatively you can use different distributions for the diagonal and off diagonal elements of $X$ in a way that will automatically satisfy the constraints): this will give you a matrix $M$ that satisfies the constraints with high probability (not surely because of the $O(\epsilon^2)$ terms).

A random unitary matrix $U$ can be generated in many ways, for example using Cayley transform: $U = (I - A)(I + A)^{-1}$ is unitary if $A$ is anti-symmetric. Therefore generate another random matrix $X$ and use its anti-symmetric part $X-X^T$ to get $U$.

You can generalize the above procedure to other functions $f$ and $g$ by expanding them in a similar way with a small parameter $\epsilon$ around $M_0$ that satisfies $f(M_0)-1=g(M_0)=0$ (because you assume $\epsilon_{1,2}$ are small). However if this condition defines a set of matrices that you want to sample from, you will have to find some parametrization that will allow you to do that (the equivalent of generating the unitary matrices in your case). There is probably no universal recipe for that so it will have to be worked out case by case.