Error using symengine in numerical solver

53 views (last 30 days)
David Keller on 24 Mar 2020 at 10:08
Commented: Walter Roberson on 28 Mar 2020 at 21:49
Hi guys
Im trying to calculate the Worst Case Expected returns according to Tüntücu and König (2004). In order to derive the robust weights, I need to solve an equation numerically. Below is the code for the calculation of the weights.
syms w_r
S = vpasolve(0==(1-((length(r_dc_s)^(-0.5))*(chi2inv(0.95,length(mu))))...
/(1*sqrt(w_r'*V*w_r)+((length(r_dc_s)^(-0.5))*(chi2inv(0.95,length(mu))))))...
*w1+w2-w_r,w_r);
w_r would be the robust weights in a vector of Nx1. When I run the code, I always get the following error message:
Error using symengine
Invalid operands.
Error in sym/privBinaryOp (line 1032)
Error in / (line 375)
X = privBinaryOp(A, B, 'symobj::mrdivide');
If I understand the error message correctly, I would have to use element-wise multiplication or division. I do not see though where this would make sense. It can be seen that only the first line of the equation would result in a 1x1 double object, so would the second. The third line contains only Nx1 vectors (w1 and w2 are the speculative and minimal weights). The only place where I can imagine it could make sense to perform a element-wise multiplication would be .*w1+w2-w_r. This leads to the same error message though..
Do you guys can help me out with this? Would be highly appreciated!
David

Walter Roberson on 28 Mar 2020 at 17:53
syms w_r
S = vpasolve(0==(1-((length(r_dc_s)^(-0.5))*(chi2inv(0.95,length(mu))))...
/(1*sqrt(w_r'*V*w_r)+((length(r_dc_s)^(-0.5))*(chi2inv(0.95,length(mu))))))...
*w1+w2-w_r,w_r);
w_r would be the robust weights in a vector of Nx1
No, when you syms and give a name, the symbolic solver always treats it as referring to a scalar. With your expression being non-scalar, you would be asking vpasolve() to find the single w_r that solves all of the elements in the vector simultaneously. vpasolve() always treats the input as simultaneous equations, and never as an array of equations to solve individually.
In the case where the desired w_r do not interact with each other, the easiest approach is:
syms w_r
S = arrayfun(@(EXPR) vpasolve(EXPR, w_r), 0==(1-((length(r_dc_s)^(-0.5))*(chi2inv(0.95,length(mu))))...
/(1*sqrt(w_r'*V*w_r)+((length(r_dc_s)^(-0.5))*(chi2inv(0.95,length(mu))))))...
*w1+w2-w_r, 'uniform', 0);
The uniform, 0 option is there because I speculate that there might not be (findable) solutions for some of the values.
In the case where the w_r values do interact with each other, then you need to create a symbolic vector:
w_r = sym('w_r', N, 1);

Show 1 older comment
Walter Roberson on 28 Mar 2020 at 20:28
Which version worked for you?
David Keller on 28 Mar 2020 at 21:46
The first option worked, the second option did not. After reading a bit about arrayfun, it is clear what I did wrong with my code at the start.
Walter Roberson on 28 Mar 2020 at 21:49
You might find that you do not need the 'uniform', 0 option. I put it in to cover the case where vpasolve() could not find a solution.

Devineni Aslesha on 27 Mar 2020 at 6:55
Hi David,
The above issue is replicated by assuming that ‘V’ is also a NX1 vector. The issue can be solved by using ‘./’ instead of ‘/’ in the vpasolve expression.
Here is a similar question for your information.

Show 1 older comment
Devineni Aslesha on 28 Mar 2020 at 6:26
Hi David
In the question, it is observed that vpasolve equation is a N*N array by considering V as an N*N matrix. This can be noticed by running the below code in which expr is an N*N array. vpasolve supports more equations than variables only for polynomial systems. For nonpolynomial equations, there is no general method of finding all solutions and vpasolve returns only one solution by default as mentioned in the doc link.
syms w_r
expr = (1-((length(r_dc_s)^(-0.5))*(chi2inv(0.95,length(mu))))...
./(1*sqrt(w_r'*V*w_r)+((length(r_dc_s)^(-0.5))*(chi2inv(0.95,length(mu))))))...
.*w1+w2-w_r;
Walter Roberson on 28 Mar 2020 at 6:39
vpasolve() simply does not support having more equations than variables for non-polynomial systems. vpasolve() uses modified newton for nonlinear cases, and that method only works for square systems (same number of equations and variables.)
David Keller on 28 Mar 2020 at 9:26
Thank you for your answers, really appreciate it. I understand the error message now, no idea though how to solve it. Seems like my approach won't work as I have hoped.
Anyways, thanks again and all the best,
David