matrixfunctions / GraphMatFun.jl

Computation graphs for matrix functions
MIT License
13 stars 0 forks source link

Fitting for graphs "linear" in some coefficients #3

Closed jarlebring closed 3 years ago

jarlebring commented 3 years ago

A graph which can be expressed as

p(x)=a+c1*f1(x)+...+ck*fk(x)

where c1,...,ck are coefficients are particularly easy to optimize. The y-coefficients in gen_general_poly_recursion (BBC-recursion) appear in that way. The other coefficients are not linear in that sense, but it can still be useful to treat these separately, since we can treat at least those coefficients perfectly.

They are easy to optimize since if we wish to fit in pointz z1,...,z2 we obtain a linear least squares problem. The setup and solution of this linear least squares problem is in the Dropbox code generators/monomial.jl function get_coeffs_linear_fit(discr,f,n0=size(discr,1);structure=:none,errtype=:abserr). It's a bad place and bad name and not good parameters. It needs clean up.

jarlebring commented 3 years ago

I see that you managed to avoid code repetition by putting crucial parts of opt_gaussnewton in separate functions. Nice work! (I deleted a redundant comment here.)

eringh commented 3 years ago

For question 5: Do we mean like minimizing |Ax-b|_1, or some type of l_1 regularization of the classical least squares problem?

Also: Do we want to try optimize it for the case of monomials? The generation of the matrix can possibly be done more efficiently, but I did not think it was worth to consider it for the initial implementation.

jarlebring commented 3 years ago

For question 5: Do we mean like minimizing |Ax-b|_1, or some type of l_1 regularization of the classical least squares problem?

I mean the situation we start with in our manuscript: max_D |p(x)-g(x)|. When p(x) has the linear structure maybe it is possible to get closer to that, rather than first relaxing to the 2-norm.

Also: Do we want to try optimize it for the case of monomials? The generation of the matrix can possibly be done more efficiently, but I did not think it was worth to consider it for the initial implementation.

Only optimize what needs to be optimized. Is this computationally demanding?

eringh commented 3 years ago

Also: Do we want to try optimize it for the case of monomials? The generation of the matrix can possibly be done more efficiently, but I did not think it was worth to consider it for the initial implementation.

Only optimize what needs to be optimized. Is this computationally demanding?

I don't think so. We are currently evaluating the graph once for the discretized points in order to create the Vandermonde-type matrix. I guess creating a proper Vandermonde matrix may be quicker, but I think we skip it then.

eringh commented 3 years ago

For question 5: Do we mean like minimizing |Ax-b|_1, or some type of l_1 regularization of the classical least squares problem?

I mean the situation we start with in our manuscript: max_D |p(x)-g(x)|. When p(x) has the linear structure maybe it is possible to get closer to that, rather than first relaxing to the 2-norm.

Ok. Not sure how to do that. I have to look into it. Perhaps something like this: https://epubs.siam.org/doi/abs/10.1137/0715015

eringh commented 3 years ago

Seems like the l1-fitting is also known as Least absolute deviations. Moreover, there are already a Julia package (LinRegOutliers) that implements it. We could also look into the JuMP package directly. Or is that too many dependencies?

eringh commented 3 years ago

No l1-optimization for now. We can open a new issue if we want to implement it later.