Hi,

I am new to Julia and trying to solve a fitting problem. I have a transition matrix that I am trying to fit. I have got some results but they are non-intutive as the probability of being in state i,i should be the greatest. How can I put such a constraint in JuMP?

The working code (not tested) is as follows:

```
function calibrate_trans_rate(filepath_pd)
# Get PD values
pd_table = zeros(Float64, 30, 7)
open(filepath_pd, "r") do fp
for row in 1:30
line = readline(fp)
pd_table[row, :] .= parse.(Float64, split(line, ","))
end
end
# Set up solver
model = Model(Ipopt.Optimizer) # need to decide what solver to use
# Set variable to be transition values
@variable(model, 0 <= x[1:8, 1:8] <= 1)
@variable(model, y[1:8, 1:8, 1:31]) # This is needed as JuMP only supports matrix operators up to quadratic terms. It doesn’t support cubic.
# Constraints on transition rates
@constraints(model, begin
default_col[i = 1:7], x[i, 8] == pd_table[1, i]
default_row[i = 1:7], x[8, i] == 0
default_state, x[8, 8] == 1
sum_of_rows[i = 1:8], sum(x[i, :]) == 1
end)
# Constraints on y {this is just a hack to get cubic terms for JuMP}
@constraints(model, begin
first_y_cons, y[:, :, 1] .== x # Or maybe set equal to x
[t = 1:30], x * y[:, :, t] .== y[:, :, t + 1]
end)
# Define statistic to minimize
@NLexpression(model, tot_sum_of_sq_diff, sum((y[k,8,time] - pd_table[1, k])^2 for k=1:7, time=1:30))
@NLobjective(model, Min, tot_sum_of_sq_diff)
JuMP.optimize!(model)
obj_value = JuMP.objective_value(model)
trans_matrix = zeros(8, 8)
for row=1:8, col=1:8
trans_matrix[row, col] = JuMP.value(x[row,col])
end
return trans_matrix
end
```

I tried adding a line in @constraints as follows:

```
curr_state[i = 1:8], x[i,i] == maximum(x[i, :])
```

This gives me the following error:

```
MethodError: no method matching isless(::JuMP.VariableRef, ::JuMP.VariableRef)
Closest candidates are:
isless(!Matched::Missing, ::Any) at missing.jl:87
isless(::Any, !Matched::Missing) at missing.jl:88
max(::JuMP.VariableRef, ::JuMP.VariableRef)@operators.jl:417
_mapreduce@reduce.jl:408[inlined]
_mapreduce_dim@reducedim.jl:318[inlined]
#mapreduce#620@reducedim.jl:310[inlined]
mapreduce@reducedim.jl:310[inlined]
_maximum@reducedim.jl:727[inlined]
_maximum@reducedim.jl:726[inlined]
#maximum#631@reducedim.jl:722[inlined]
maximum@reducedim.jl:722[inlined]
macro expansion@rewrite.jl:227[inlined]
macro expansion@macros.jl:440[inlined]
(::Main.workspace192.var"#7#15"{JuMP.Model,Array{JuMP.VariableRef,2}})(::Int64)@macro.jl:183
#26@container.jl:70[inlined]
iterate@generator.jl:47[inlined]
collect(::Base.Generator{JuMP.Containers.VectorizedProductIterator{Tuple{Base.OneTo{Int64}}},JuMP.Containers.var"#26#27"{Main.workspace192.var"#7#15"{JuMP.Model,Array{JuMP.VariableRef,2}}}})@array.jl:686
map@abstractarray.jl:2188[inlined]
container@container.jl:70[inlined]
container@container.jl:65[inlined]
macro expansion@macros.jl:91[inlined]
calibrate_trans_rate(::String)@Other: 25
top-level scope@Local: 1
```

Any ideas?

The problem that I am trying to solve is here: