""" Calculates the Tikhonov regularized estimate min. ||Ax-b|| + λ||x|| using SVD for a vector with different values of λ `x, normE, normX = dsvd(A,U,S,V,b,λ)` `x, normE, normX = dsvd(A,b,λ)` """ function dsvd(A,U,S,V,b,λ) if minimum(λ) < 0 error("Illegal regularization parameter λ") end UTb = U'b N = size(U,1) n = size(V,1) nl = length(λ) x = zeros(n,nl) normE = zeros(nl) normX = zeros(nl) print("dsvd: ") for (i,λi) in enumerate(λ) x[:,i] = V*diagm(S./(S.^2+λi))*UTb normX[i] = norm(x[:,i]) normE[i] = norm(A*x[:,i] - b) print("[\$(round(i/length(λ),2)),\$(round(λi,2))] ") end println("Done ") return x, normE, normX end function dsvd(A,b,λ) if minimum(λ) < 0 error("Illegal regularization parameter λ") end U,S,V = svd(A) dsvd(A,U,S,V,b,λ) end """ Lcurve(normE, normX, λ) Plots the L-curve. normE and normX are obtained from, e.g., `dsvd` """ function Lcurve(normE, normX, λ) plot(normE,normX,xscale=:log10,yscale=:log10,m=:o, xlabel="RMSE", ylabel="||k||", title="L-curve") annotations = [(normE[i],normX[i],"λ=\$(round(λ[i],8))") for i in 1:length(λ)] annotate!(annotations) end