SciML / DataDrivenDiffEq.jl

Data driven modeling and automated discovery of dynamical systems for the SciML Scientific Machine Learning organization
https://docs.sciml.ai/DataDrivenDiffEq/stable/
MIT License
404 stars 57 forks source link

optimal thresholding for truncated SVD ? #315

Open yewalenikhil65 opened 2 years ago

yewalenikhil65 commented 2 years ago

Should we be considering optimal thresholding in truncated SVD (maybe as another dispatch for truncated_svd) ? i think some DMD algorithms with noisy data try to employ these tricks. Paper The Optimal Hard Threshold for Singular Values is 4/sqrt(3) is a useful reference

AlCap23 commented 2 years ago

Ah, nice! I did not know about this ( albeit it makes total sense here. )

The algorithm is already implemented here, since it is used for denoising of the data in the sparse regression.

It is definitely possible, but I would need to think about the dispatch.

yewalenikhil65 commented 2 years ago

Ah, nice! I did not know about this ( albeit it makes total sense here. )

The algorithm is already implemented here, since it is used for denoising of the data in the sparse regression.

It is definitely possible, but I would need to think about the dispatch.

Ohh.. I didn't knew that it's already there in DataDrivenDiffEq.jl. I will try to see how we can add a dispatch and if it works as expected

yewalenikhil65 commented 2 years ago

@AlCap23 , we can try to do something similar to usage of fbDMD with optimal_svht https://github.com/amir-cardiolab/DMD_DA_blood_flow/blob/53f1088a9dd003cf0d2d5c6c34d2c0d0dfcf1cfd/ROM-KF/ROM_KF.m#L57

This is the code i followed for FBDMD PR as well, and is based on paper 2 mentioned in the PR of FBDMD

My suggested dispatch as follows

# truncated_svd ,based on  optimum SVHT 
function truncated_svd(A::AbstractMatrix{T}) where T <: Number
    U, S, V = svd(A)  # This is already equal to "econ" from MATLAB
    truncation = optimal_svht(size(S,1), size(S,1))*median(S)
    r = length(findall(x -> x > truncation, S));
    U = @view U[:, 1:r]
    S = @view S[1:r]
    V = @view V[:, 1:r]
    return U, S, V
end

Another thing to think over should be how should such dispatch be facilitated for existing DMD functions in the package? My suggestion is to add a kwarg optimal = false for existing DMDSVD n DMDINV functions , which by default should call existing truncated_svd, and if specified optimal= true should call new dispatch of truncated_svd based on optimal_svht function written above!