mm_roots enlarging argument limitation

[Eileen Witowski]

I am using the mm_root function and I noticed that the the the number of arguments is limited by 10. However I need to use more arguments. So I looked up the source Code of mm_root and tried to change the code for my needs. I want to pass 75 arguments to f, not only 10. So I changed the relevant line in the code from

local opts "o1, o2, o3, o4, o5, o6, o7, o8, o9, o10"

to:
local opts "o1, o2, o3, o4, o5, o6, o7, o8, o9, o10, o11, o12, o13, o14, o15, o16, o17, o18, o19, o20, o21, o22, o23, o24, o25, o26, o27, o28, o29, o30, o31, o32, o33, o34, o35, o36, o37, o38, o39, o40, o41, o42, o43, o44, o45, o46, o47, o48, o49, o50, o51, o52, o53, o54, o55, o56, o57, o58, o59, o60, o61, o62, o63, o64, o65, o66, o67, o68, o69, o70, o71, o72, o73, o74, o75"

Before changing I got the following error:
mm_callf_setup(): 3001 expected 4 to 16 arguments but received 81

After changing the Code, I now get another error, which I do not understand:
mm_callf_setup(): 3001 expected 2 to 12 arguments but received 78

Trying to understand that error I made stepwise changes, e.g. I added 10 additional locals (

local opts "o1, o2, o3, o4, o5, o6, o7, o8, o9, o10, o11, o12, o13, o14, o15, o16, o17, o18, o19, o20") and so on, up to

local opts "o1, o2, o3, o4, o5, o6, o7, o8, o9, o10"

to:

local opts "o1, o2, o3, o4, o5, o6, o7, o8, o9, o10, o11, o12, o13, o14, o15, o16, o17, o18, o19, o20, o21, o22, o23, o24, o25, o26, o27, o28, o29, o30, o31, o32, o33, o34, o35, o36, o37, o38, o39, o40, o41, o42, o43, o44, o45, o46, o47, o48, o49, o50, o51, o52, o53, o54, o55, o56, o57, o58, o59, o60, o61, o62, o63, o64, o65, o66, o67, o68, o69, o70, o71, o72, o73, o74". The Programm returned the error I expected: mm_callf_setup(): 3001 expected 4 to 80 arguments but received 81

That's why I don't understand the error after allowing 75 locals.

I will be very happy about useful hints.

[William Lisowski]

Welcome to Statalist.

I'll start by saying I'm not the Mata expert around here, but I do have a possible workaround similar to what I've done elsewhere.

Assuming that your arguments are scalars, have you considered rewriting the function you are finding the root of, to which you need to pass 75 arguments, to instead accept a single vector containing those 75 arguments?

mata
mata clear

function myfunc(x, a) return(x^2 - a)
a = 2/3
a
mm_root(x=., &myfunc(), 0, 1, 0, 1000, a)
x

function myfund(x, b) return(x^2 - b[1,1])
b = J(2,1,.)
b[1,1] = 2/3
b
mm_root(x=., &myfund(), 0, 1, 0, 1000, b)
x

end

: mata clear

: 
: function myfunc(x, a) return(x^2 - a)

: a = 2/3

: a
  .6666666667

: mm_root(x=., &myfunc(), 0, 1, 0, 1000, a)
  0

: x
  .8164965809

: 
: function myfund(x, b) return(x^2 - b[1,1])

: b = J(2,1,.)

: b[1,1] = 2/3

: b
                 1
    +---------------+
  1 |  .6666666667  |
  2 |            .  |
    +---------------+

: mm_root(x=., &myfund(), 0, 1, 0, 1000, b)
  0

: x
  .8164965809

: 
: end

[Eileen Witowski]

Thanks for your answer! I think this will probably not work for my function.
It would have been smart to post my function right away. I will give a short example, but notice I would like to use a1,...,a37, b1,...,b37, V.

b = 608;
z = 208;
roots = J(b,z,0)

function myfunc(root, a1, a2, a3, a4, a5, a6, b1, b2, b3, b4, b5, b6, V) return( a1/((1+root)^b1) + a2/((1+root)^b2) + a3/((1+root)^b3) + a4/((1+root)^b4) + a5/((1+root)^b5) + a6/((1+root)^b6) - V )

mm_root(root=., &myfunc(), 0, 1, 0, 1000, c1[1], c2[1], c3[1], c4[1], c5[1], c6[1], t1[1], t2[1], t3[1], t4[1], t5[1], t6[1], P[1,1])
roots[1,1] = root

where P is a Matrix of size 208x608 and c as well as t are vectors of length 608.

[William Lisowski]

Either of the following approaches replaces your 13 scalar arguments with 3 arguments: 2 vectors of 6 scalars and one further scalar.

function myfunc(root, a, b, V) {
return( a[1]/((1+root)^b[1]) + a[2]/((1+root)^b[2]) + a[3]/((1+root)^b[3]) + a[4]/((1+root)^b[4]) + a[5]/((1+root)^b[5]) + a[6]/((1+root)^b[6]) - V )
}

function myfunc(root, a, b, V) {
a1 = a[1]
a2 = a[2]
a3 = a[3]
a4 = a[4]
a5 = a[5]
a6 = a[6]
b1 = b[1]
a2 = b[2]
a3 = b[3]
a4 = b[4]
a5 = b[5]
a6 = b[6]
return( a1/((1+root)^b1) + a2/((1+root)^b2) + a3/((1+root)^b3) + a4/((1+root)^b4) + a5/((1+root)^b5) + a6/((1+root)^b6) - V )
}

Then you create your vector arguments and call mm_root.

avec = (c1[1]\c2[1]\c3[1]\c4[1]\c5[1]\c6[1])
bvec = (t1[1]\t2[1]\t3[1]\t4[1]\t5[1]\t6[1])
mm_root(root=., &myfunc(), 0, 1, 0, 1000, avec, bvec, P[1,1])

Note that this code is untested and is intended only to indicate an approach that may work for you, not to be a complete solution.

[Eileen Witowski]

Even though I'm not getting results that seem to be right yet, the approach seems to work just fine. Thank you!

[Eileen Witowski]

The mm_root function

issues return code rc with rc!=0 indicating that no valid solution has been found. Since I have a big data set and therefore a lot of rc's, I would like to select those where rc!=0 and then analyze the underlying problem…
Did anybody make experiences? I would be glad About hints to realize that.
Thanks!

[William Lisowski]

Without some idea of the code you are using to, apparently, run mm_root repeatedly - on different subsets of your data? on the same data with different sets of parameter values? - it's difficult to recommend an appropriate method from among the many that will accomplish what you describe.

It's probably time to remind you of the Statalist FAQ linked to from the top of the page, as well as from the Advice on Posting link on the page you used to create your post. If you haven't done so already, you should take the time to review it, noting especially sections 9-12 on how to best pose your question.

[Eileen Witowski]

Let me specify my question. I defined 37 functions that look just like this:
function myfunc37(root, a, d, b, ip, V) {
b1 = b[1]
b2 = b[2]
b3 = b[3]
b4 = b[4]
b5 = b[5]
b6 = b[6]
b7 = b[7]
b8 = b[8]
b9 = b[9]
b10 = b[10]
b11 = b[11]
b12 = b[12]
b13 = b[13]
b14 = b[14]
b15 = b[15]
b16 = b[16]
b17 = b[17]
b18 = b[18]
b19 = b[19]
b20 = b[20]
b21 = b[21]
b22 = b[22]
b23 = b[23]
b24 = b[24]
b25 = b[25]
b26 = b[26]
b27 = b[27]
b28 = b[28]
b29 = b[29]
b30 = b[30]
b31 = b[31]
b32 = b[32]
b33 = b[33]
b34 = b[34]
b35 = b[35]
b36 = b[36]
b37 = b[37]
return( a/((d+root)^b1) + a/((d+root)^b2) + a/((d+root)^b3) + a/((d+root)^b4) + a/((d+root)^b5) + a/((d+root)^b6) + a/((d+root)^b7) + a/((d+root)^b8) + a/((d+root)^b9) + a/((d+root)^b10) + a/((d+root)^b11) + a/((d+root)^b12) + a/((d+root)^b13) + a/((d+root)^b14) + a/((d+root)^b15) + a/((d+root)^b16) + a/((d+root)^b17) + a/((d+root)^b18) + a/((d+root)^b19) + a/((d+root)^b20) + a/((d+root)^b21) + a/((d+root)^b22) + a/((d+root)^b23) + a/((d+root)^b24) + a/((d+root)^b25) + a/((d+root)^b26) + a/((d+root)^b27) + a/((d+root)^b28) + a/((d+root)^b29) + a/((d+root)^b30) + a/((d+root)^b31) + a/((d+root)^b32) + a/((d+root)^b33) + a/((d+root)^b34) + a/((d+root)^b35) + a/((d+root)^b36) + (a+ip)/((d+root)^b37) - V )
}

After defining my functions I am using mm_root to solve for the variable "root". I'm using a case-by-case analysis. To make it short I'm just giving you parts of the code here as well:
i = 0
for(k=1; k<209; k++){
for(j=1; j<609; j++){
i = j + (k-1)*608
if (cf[i]!=. & t1[i]!=. & t2[i]!=. & t3[i]!=. & t4[i]!=. & t5[i]!=. & t6[i]!=. & t7[i]!=. & t8[i]!=. & t9[i]!=. & t10[i]!=. & t11[i]!=. & t12[i]!=. & t13[i]!=. & t14[i]!=. & t15[i]!=. & t16[i]!=. & t17[i]!=. & t18[i]!=. & t19[i]!=. & t20[i]!=. & t21[i]!=. & t22[i]!=. & t23[i]!=. & t24[i]!=. & t25[i]!=. & t26[i]!=. & t27[i]!=. & t28[i]!=. & t29[i]!=. & t30[i]!=. & t31[i]!=. & t32[i]!=. & t33[i]!=. & t34[i]!=. & t35[i]!=. & t36[i]!=. & t37[i]!=.) {
bvec = (t1[i]\t2[i]\t3[i]\t4[i]\t5[i]\t6[i]\t7[i]\t8[i]\t9[i]\t10[i]\t11[i]\t12[i]\t13[i]\t14[i]\t15[i]\t16[i]\t17[i]\t18[i]\t19[i]\t20[i]\t21[i]\t22[i]\t23[i]\t24[i]\t25[i]\t26[i]\t27[i]\t28[i]\t29[i]\t30[i]\t31[i]\t32[i]\t33[i]\t34[i]\t35[i]\t36[i]\t37[i])
mm_root(root=., &myfunc37(), -1, 100, 0.001, 10000, cf[i], dr[i], bvec, nw[i], P[k,j])
roots[k,j] = root
}
…….
……..
}
}

Note roots is a 208x608 matrix.
Running the programm I noticed some rc!=0. Now I would like to locate and analyze those.

[William Lisowski]

Consider the following, where I chose an argument that causes mm_root to have a non-zero return code.

mata
mata clear
function myfund(x, b) return(x^2 - b[1,1])
b = J(2,1,.)
b[1,1] = -5
rc=mm_root(x=., &myfund(), 0, 1, 0, 1000, b)
rc
x
end

: mata clear

: function myfund(x, b) return(x^2 - b[1,1])

: b = J(2,1,.)

: b[1,1] = -5

: rc=mm_root(x=., &myfund(), 0, 1, 0, 1000, b)

: rc
  2

: x
  0

: end

By having replaced

mm_root(x=., &myfund(), 0, 1, 0, 1000, b)

with

rc=mm_root(x=., &myfund(), 0, 1, 0, 1000, b)

the return code is no longer displayed, but instead is stored in rc, just as the root is stored in x. So you can create another 208x608 matrix to store the return code, in the same way that you store the root. And you can afterwards find the nonzero elements in the matrix of return codes.