Hi,
I’m trying to run a nonlinear regression fit for a perturbation function, but it seems to only fit one of my two parameters. The perturbation function is for modeling grating efficacy curves, with wavelength as the x-axis and efficacy as the y-axis. This function utilizes a library of grating models I’ve created where each model is set at specific a Duty Cycle value and spans several Depth values. I’m trying to use LsqFit to best fit my two parameters, Depth and Duty Cycle (D and DC), with the wave index set as the one-element predictor variable.
Below are my three functions, one that creates the perturbation function from my library of grating models, one that reads in an “unknown” efficiency curve, and one that runs the nonlinear regression. I’ve simplified the functions here by using a sample of models covering the Depth range from 1 to 2 um and wavelengths 400 to 740 nm. I’m also using fewer base models to create the perturbation function (DC_list in code below) and I saved the base model information used to create the perturbation as “.jl” files (these are actually “.jld” files that I renamed to “.jl” file so I could upload them at the very bottom here, they run as “.jld” files using the load() function). These files include a list of depths, a list of wavelengths, and several matrices of efficiencies for specific Duty Cycle setups. The efficiency matrices are sliced depending on the desired Depth that then gives us the final efficiency curve. I use wave index rather than just wavelength in the nonlinear regression since later I will need to look at certain wavelength ranges for real data.
My problem is this: When I run my Lsqfit nonlinear regression function to find the best fitting Depth and Duty Cycle values, the Duty Cycle parameter will vary, but the Depth value will stay equal to my initial parameter guess. In other words, when I run curve_fit with the initial parameter value p0 = [D_0, DC_0], the best fit will always have the Depth parameter = D_0, while the Duty Cycle parameter will vary, leading to incorrect parameter fitting (example at bottom after code).
To test the function, I’m trying to fit one of the base models use to create the perturbation model, since this should give me the exact/very close to the true Depth and Duty Cycle values. I’ve set up the base model range (DC_list) right now to fit DC = 28.0 (this is matrix file “z_matrix_DC28.0.jl”), with the D parameter able to be anything between 1.2 to 1.6. I have tested the nonlinear regression function for several different parameter values using my full library of grating models, but I get the same result where the best fit Duty Cycle parameter will differ from DC_0, but the best fit Depth parameter will always equal the initial parameter D_0.
Here is the perturbation function:
using CubicSplines
using LsqFit
using JLD
function Pert_function(mod_wave_index, mod_D, mod_DC)
# Create Depth and Duty Cycle list used to create the perturbation fit
DC_list = [23.0, 26.5, 27.0, 27.5, 28.0, 28.5, 29.0, 29.5, 50.0]
num_list = collect(1:1:length(DC_list))
D_list = collect(1.1:0.01:1.7)
path_to_all_mods = "C:\\Users\\norma\\Desktop\\G_MODS\\Model_x_d_and_z_data\\"
cd(path_to_all_mods)
z_list = []
for n in num_list
#Bring in Base Model data
d = load("d_vector_ALL.jl")["data"]
z = load("z_matrix_DC"*string(DC_list[n])*".jl")["data"]
# Slice matrix for specific value of depth (mod_D)
dδ = abs.(d .- mod_D) #<------------------- 1st parameter for the nonlinear regression fitting
idδ = findall(minimum(dδ) .== dδ)[1]
z_good = z[:,idδ]
append!(z_list, [z_good])
end
####### Spline Perturbation Fitting
f_all_list = []
for i in mod_wave_index
effi_list = []
for n in num_list
z_list_add = z_list[n][i]
append!(effi_list, z_list_add)
end
effi_list = Float64.(effi_list)
#Create spline fit for perturbation function at every wavelength
f = CubicSpline(DC_list, effi_list)
#Final model
f_v = f(mod_DC)#<------------------- 2nd parameter for the nonlinear regression fitting
# Save perturbation fit efficancies
append!(f_all_list, f_v)
end
# Change to Float64
f_all_list = Float64.(f_all_list)
return f_all_list
end
Here is the function that reads in the “unknown” data:
function unknown_data(UK_wave_index, UK_D, UK_DC)
# Bring in "unknown" data and prep for comparisons
path_to_UK_to_fit = "C:\\Users\\norma\\Desktop\\G_MODS\\Model_x_d_and_z_data\\"
cd(path_to_UK_to_fit)
x = load("x_vector_ALL.jl")["data"]
d = load("d_vector_ALL.jl")["data"]
z = load("z_matrix_DC"*string(UK_DC)*".jl")["data"]#<--------------- 2nd parameter set for "unknown" data set
dδ = abs.(d .- UK_D) #<--------------- 1st parameter set for "unknown" data set
idδ = findall(minimum(dδ) .== dδ)[1]
z_list = z[:,idδ]
wave_all_list = []
effi_list = []
#Pick wavelength ranges/values
for i in UK_wave_index
wave_used = x[i]
eff_list = z_list[i]
append!(effi_list, eff_list)
append!(wave_all_list, wave_used)
end
# Change to Float64
effi_list = Float64.(effi_list)
wave_all_list = Float64.(wave_all_list)
return wave_all_list, effi_list
end
Finally, here is the nonlinear regression function that should fit the “unknown” data to the perturbation function:
function nonlinear_reg_fitting_UK_model(UK_wave_index, UK_D, UK_DC, d_dc_guess)
#Set UK_wave_index = collect(1:341) for full wave range
UK_wave_values, UK_effi_values = unknown_data(UK_wave_index, UK_D, UK_DC)
m(mod_wave_index, d_dc) = Pert_function(mod_wave_index, d_dc[1], d_dc[2]) # d_dc = [D, DC]
fit = curve_fit(m, UK_wave_index, UK_effi_values, d_dc_guess, show_trace=true)
print(" Real D and DC: ")
print([UK_D, UK_DC])
print(" Guess D and DC: ")
print(d_dc_guess)
print(" Best fit D and DC: ")
print(fit.param)
end
As an example, if I input:
nonlinear_reg_fitting_UK_model(collect(1:341), 1.4, 28.0, [1.2, 27.0])
into julia, this will return:
0 2.833372e+00 NaN
* lambda: 10.0
1 2.459562e+00 4.004799e-01
* g(x): 0.40047990579591175
* lambda: 1.0
* dx: [-0.0, 0.4295621750066164]
2 1.127579e+00 2.367417e-01
* g(x): 0.23674172253530493
* lambda: 0.1
* dx: [-0.0, 2.117552181743514]
3 3.510501e-01 5.600239e-02
* g(x): 0.05600239492178097
* lambda: 0.010000000000000002
* dx: [-0.0, 3.1952790213651037]
4 2.594502e-01 4.757896e-03
* g(x): 0.004757896119312277
* lambda: 0.0010000000000000002
* dx: [-0.0, 1.6641849999546636]
5 2.557867e-01 1.260775e-04
* g(x): 0.00012607754808589974
* lambda: 0.00010000000000000003
* dx: [-0.0, 0.3864222083774912]
6 2.557155e-01 2.049823e-06
* g(x): 2.04982345044857e-6
* lambda: 1.0000000000000004e-5
* dx: [-0.0, 0.05514760123727751]
7 2.557143e-01 3.370836e-08
* g(x): 3.370835833004503e-8
* lambda: 1.0000000000000004e-6
* dx: [-0.0, 0.007235153001275311]
8 2.557143e-01 5.605855e-10
* g(x): 5.605854789887964e-10
* lambda: 1.0000000000000005e-7
* dx: [-0.0, 0.0009398515211527311]
9 2.557143e-01 9.356143e-12
* g(x): 9.35614322800556e-12
* lambda: 1.0000000000000005e-8
* dx: [-0.0, 0.00012193581672843849]
10 2.557143e-01 1.565193e-13
* g(x): 1.5651928082022733e-13
* lambda: 1.0000000000000005e-9
* dx: [-0.0, 1.581737817222516e-5]
Real D and DC: [1.4, 28.0] Guess D and DC: [1.2, 27.0] Best fit D and DC: [1.2, 34.85646094540199]
Where we can see that Lsqfit is not varying the D parameter and therefor best fit for Depth equals the initial guess value of 1.2. If I set the Depth initial parameter as correct value of 1.4 for this case (d_dc_guess = [1.4, 27.0]) we get a DC best fit of about 28.0, meaning that the nonlinear regression is working for the DC parameter, just not for the D parameter.
I can not figure out what the problem here is, and I will take any suggestions people may have. I just need both parameters to be changing, not just the one. Is Lsqfit the best way to go about the nonlinear regression for my perturbation function? Am I missing something completely in my final function? Thanks for any help in advance!
DATA USED:
Here are the files I used that were located in my Model_x_d_and_z_data folder.
d_vector_ALL.jl (5.1 KB)
x_vector_ALL.jl (7.0 KB)
z_matrix_DC23.0.jl (273.4 KB)
z_matrix_DC26.5.jl (273.4 KB)
z_matrix_DC27.0.jl (273.4 KB)
z_matrix_DC27.5.jl (273.4 KB)
z_matrix_DC28.0.jl (273.4 KB)
z_matrix_DC28.5.jl (273.4 KB)
z_matrix_DC29.0.jl (273.4 KB)
z_matrix_DC29.5.jl (273.4 KB)
z_matrix_DC50.0.jl (273.4 KB)