kr-colab/product-space-fem

product-space-fem's Issues

gradient verification

We need a method that compares the adjoint gradient dJ/dm to the difference quotient of the loss functional. By Taylor's theorem the convergence should be
|J(m + h m_) - J(m) - h dJ/dm • m_ | -> 0 at O(h^2)
where m_ is an arbitrary function in the same space as m.

name!

conda environment works with flexible but not strict channel priority

I'd read that setting channel_priority to strict would make things faster/better; however, that doesn't work with our environment here. Maybe that's not a problem, but I'm recording the fact in case it helps clear up later errors.

I think that the thing that made it fail was channel priority; but I've still got to verify that.

For the record, here was the error:
err.txt

get contours

using gdal-contour

outline easy-mode workflow

Input data:

• data frame of spatial locations of samples and labels
• data frame of heterozygosities and diversities
• data frame of coordinates of corners of range polygon, counter-clockwise

Input parameters:

• epsilon
• mesh resolution
• regularization strengths (four of these?)
• filenames of the three input data files
• output directory basename

Steps:

0. (make the range polygon)
1. run the method
2. (do crossvalidation?)
3. tweak plots (?)

Output:

1. optimization diagnostics: history of parameter values and loss values, messages
2. optimal ellipse and vector fields:
   • images (png or pdf)
   • optimized equation (including mesh, function space, optimal values, etc) (if pickle-able??)
3. crossvalidation error (?)
   • plot
   • values in csv

Advanced parameters:

• scipy-minimize tolerances and other options

To do crossvalidation:

• command to make one crossvalidation subdirectory
• script/notebook to make bash script to run lots of crossvalidation things
• script/notebook to collate results from many subdirectories

How to get the range polygon?

• raster file -> gdal-contour -> polygon
• ellipse/rectangle etc with meshr -> ??
• should standardize coordinates, at least internally?
• "convex hull" can get counter-clockwise points

plotting ProductFunction

We need a plot() method in the ProductFunction class with a slice as an argument. For lower dimensions slice can be None and plot the full (e.g. 2D) function. Alternatively we can plot a slice of the function by providing a point p in one of the marginal spaces, which would plot F(x) = f(p, x).

scipy >= 1.8.1 causes a segfault

@silastittes has tracked it down to a change that happened between 1.8.0 and 1.8.1. This is puzzling, since not much changed there, except for their unuran submodule. Debugging this seems daunting, but it might be a legit scipy bug.

Here is the script that causes the segfault: https://github.com/kr-colab/product-space-FEM/blob/daf7a4580dba504b622ca240c76d3c9949449f14/debug.py

`InverseProblem.optimize( )` both modifies its input in place and returns it

This is confusing since it's not obvious that it's updated it in place; the easiest solution would be to just not return it at all.

Ellipse field plot

Right now any Function with a value dimension under 3 has a .plot() method.
E.g. the covariance matrix in the spatial coalescence equation does not have a `.plot()' method.

As a special case, the unconstrained parameterization of the covariance matrix should have a plot showing ellipses at each point. Point-wise, the ellipses are shaped by the eigenvectors of the covariance matrix.

ensure the `pf.SpatialData` class gets data in the right order

Right now the SpatialData class takes data and points, assuming that the data are in upper-triangular order (or, maybe lower-triangular?). This should be more robust; possibly we can just move over and merge the SpatialDivergenceData class I've defined in inference/?