I would like to optimize a set of 3 variable vectors to minimize the total variance in a dataset. I think my problem is in passing variable and parameter vectors into my user-defined objective.

I found this question on StackOverflow, but I don’t know how to use it with several variable and parameter vectors.

Here’s the commented code with sample data, maybe self-explanatory?

```
# User-defined objective:
# 👇 This is where I don't know what I'm doing 👇
# I don't know how to pass the variables and parameters,
# and I don't know how to access them from within the function,
# and I need γ to have an integer constraint, I plan to simply use floor() once I can access the variable within this objective.
function xy_objective(α..., β..., γ..., Y_0..., nan_pad_fraction)
# 👆 This is where I don't know what I'm doing 👆
n_p, n_r, n_t = size(Y_0)
n_nans = convert(Int,round(nan_pad_fraction*n_p))
# 1. circshift each column by its γ (by same amount for all targets)
nan_pad = replace!(zeros(n_nans,n_r,n_t), 0=>NaN)
Y_0_padded = vcat(nan_pad, Y_0, nan_pad)
Y_0_shift = Y_0
for i_r = 1:n_r
Y_0_shift[:,i_r,:] = circshift(Y_0_padded[:,i_r,:], γ[i_r])
end
# 2. take nanmean of each target to reduce variance against target mean
ȳ = reshape(nanmean(Y_0_shift,dims=2),n_p + 2*n_nans, n_t)
# 3. calculate cost
χ²_t = zeros(n_t)
for i_t = 1:n_t
χ²_t =
sum(Y_0_shift[:,:,i_t] -
(reshape(α,(1,n_r)) .+ ȳ[:,i_t]*reshape(β,(1,n_r))).^2)
end
χ² = sum(χ²_t)
return χ²
end
#################################################################################
# Optimization:
# Sample data
n_pts = 1000 # number of data sample points per experiment
n_reps = 10 # number of experimental replicates
n_targ = 2 # number of experimental targets (simultaneously acquired in each rep)
Y_0 = rand(n_pts, n_reps, n_targ) # data matrix
# Goal: reduce total variance by:
# adding α to each column
# scaling scaling by β for each column
# shifting by γ for each target
# Pad edges with NaN values so data can shift
np_frac = 0.1
# e.g. a vector with 100 elements along the first dimension will have 100*np_frac NaNs added to beginning and end to become 100*(1+np_frac) elements long in the first dimension
# Build model
m = Model(Ipopt.Optimizer)
register(m,
:xy_objective,
n_reps*(2*n_targ + 1),
xy_objective; autodiff=true)
α_test, β_test, γ_test = nothing, nothing, nothing
nan_pad_fraction, Y_0_prob = nothing, nothing
# α: additive variable, 1 for each replicate for each target
@variable(m,
α_test[i=1:n_reps,j=1:n_targ] == 0)
# β: scaling variable, 1 for each replicate for each target
@variable(m,
β_test[i=1:n_reps,j=1:n_targ] == 1)
# γ: shifting variable, 1 for each replicate (shift all targets by same amount)
@variable(m,
-convert(Int,floor(np_frac*n_pts)) <=
γ_test[i=1:n_reps] <=
convert(Int,floor(np_frac*n_pts)) ) # also need integer constraint
# Fraction of data length to allow for shifting
@NLparameter(m,
nan_pad_fraction == np_frac)
# Starting point for experimental data
@NLparameter(m,
Y_0_prob[i=1:n_pts,j=1:n_reps,k=1:n_targ] == Y_0[i,j,k])
# 👇 This is where I don't know what I'm doing 👇
# I don't know how to pass the variables and parameters,
# and I don't know how to access them from within the function.
@NLobjective(m,
Min,
xy_objective(α_test..., β_test..., γ_test..., Y_0_prob..., nan_pad_fraction))
# 👆 This is where I don't know what I'm doing 👆
optimize!(m)
```