I have to estimate the coefficients $c_1,...,c_n$ of the relationship $$ Y = c_1X_1+c_2X_2+\cdots + c_nX_n $$ So the distribution of $Y$ is given through a linear combination of other random variables. Further there is the fact that for any sample of $X_1,...,X_n$ it holds that $$ 0\leq\sum_{i=1}^nX_i<1 $$ I have ways to generate samples for $Y,X_1,...,X_n$, giving me a linear system of equations. One way to estimate the coefficients $c_1,...,c_n$ in that case would be through least squares. But I'd like to do better than that.
Another alternative would be to use moment matching (with centralized moments), e.g., $$ E[(Y-E[Y])^m] = E[(c_1X_1+c_2X_2+\cdots + c_nX_n-E[c_1X_1+c_2X_2+\cdots + c_nX_n])^m] $$ which would add further equations for higher moments with $m>1$.
This however leads to a system of equations which is non-linear in the parameters $c_1,...,c_n$.
Are there any principled ways to solve such non-linear systems?
The best idea I had so far (for $m=1,...,N$) is to find a numerical solution through optimizing the squared loss of the LHS and the RHS of the equation through gradient descent.
Do you have any better ideas or principled approaches to suggest?