Problems with the parametric Recurrence Table

I am trying to construct a minimization program where I have to minimize a function with respect to a matrix. This matrix enters into the optimization program in a convoluted manner together with data. My problem is that I cannot construct this part of the function.

What I have tried:

First generate some fake data:

A = {{5, 4, 1}, {4, 7, 3}, {1, 3, 5}};
B = {{3, 2, -1}, {2, 5, 1}, {-1, 1, 3}};
Xsam = RandomVariate[MultinormalDistribution[{0, 0, 0}, B], 1000];
Ysam = RandomVariate[MultinormalDistribution[{0, 0, 0}, A], 1000];

Generate the parametric matrix

M0 = Array[Subscript[a, ##] &, {3, 3}];

The first step in creating the final function of M, Xsam, Ysam is:

r[X_,Y_,M_]:=Exp[X.M.Transpose[Y]];

Then, for

pX=1/1000;
pY=1/1000;
r1=r[Xsam,Ysam,M0]

I would like to take the 10th of iteration of a recurrence relationship defined in this table (I already use "Take") which should be a function of M0

Take[RecurrenceTable[{b[n + 1] == pY/(Transpose[r1].a[n]),
a[n + 1] == pX/(r1.b[n]), a == Exp[-ConstantArray[1, Length [Xsam]]],
b == Exp[-ConstantArray[1, Length [Xsam]]]}, {a, b}, {n, 1,
10}], {10}]

However it does not work (I get the message "RecurrenceTable: The function appears with no arguments").

If I use instead of M0 some numerically valued matrix, it does work (even if surprisingly slow, even if on 32G RAM PC).

Thank you for any hints as to how to get to the 10th iteration (and make it defined as a function of elements of M0) on which I will be later optimizing on.

• Please explain in more detail what you actually want to do. This does not make sense. Do you really want to carry the symbolic entries along through all computations? Then I am tempted to say that won't work because the complexity of the symbolic expression will blow up pretty fast. – Henrik Schumacher Aug 18 at 3:59
• Thank you very much for your comment. I essentially need to fit an equation in M0 (matrix) subject to two other equation constraints (which depend on M0 and also some functions b[x] and a[y] which are themselves unknown. So I started from iteratively expressing b(x) and a(y) through M0 and wanted to run NMinimize to find M0 that minimize my fitting equation. But I realize that it won't fly, indeed it would be too complex to keep track of M0 as a symbolic expression. I will try to fit directly all 3 pieces (M0, b(.) and a(.)) through a different code.... – Kass Aug 18 at 17:20