# How should I solve the following simultaneous equations?

I have the set of simultaneous equations below from the paper Parameter estimation for 3-parameter generalized pareto distribution by the principle of maximum entropy (POME) by VP Singh and H Guo:

$$\sum_{i=1}^{n} \frac{(x_i - c)/\hat{b}}{1-\hat{a}(x_i-c)/\hat{b}} = \frac{n}{1-\hat{a}}$$

$$\sum_{i=1}^{n} ln[1-\hat{a}(x_i-c)/\hat{b}] = -n\hat{a}$$

In this problem $c$ is set by me. $x$ is a vector of returns $\geq c$ and $n$ is the number of elements in $x$. Does anyone have any recommendations for how to solve for $a$ and $b$? I generally use bisection as my go to root solver but I have never used it to find two variables.

For what it's worth I am programming in C++ per my advisor's request but that is the next question.

• On an answer to an earlier question of yours here I gave the quite explicit advice to minimize the negative log-likelihood by directly solving the minimization problem using standard tools, rather than taking the step in the paper you point to of taking derivatives and then trying to make it a root-finding problem (which they give no advice on). Since you're fairly likely to get the same advice again, you should explain why you don't do it. Nov 3, 2015 at 21:43
• I spent the weekend working on the minimization problem you reference however I couldn't find standard tools that both worked well with visual studio and minimization problems with a vector in them. Nov 3, 2015 at 21:50

On thinking about this for a bit, I think I see how to solve this problem directly without any need for iteration.

I'll refer to the Wikipedia pages for the Generalized Pareto and the Pareto. You'll now need to navigate between three sets of notation, but the critical one is $a=\xi=1/\alpha$. Note also that in the discussion below $x_{(1)}$ means the first order statistic (the smallest observation).

(Note added much later: here I seem to be assuming $a>0$.)

1. Given $c$ and assuming we have some suitable $\hat{b}$, let $z_i = \frac{x_i-c}{\hat{b}}$ (you'll need these to be such that all the $x_i>c$ and all the $z_i>1$). You now have a standard generalized Pareto, which has only one parameter:

$f_a(z) = (az+1)^{-\frac{a+1}{a}}$

You could of course try univariate optimization on this new problem, but there's a direct solution.

Use the connection between the Pareto and the Generalized Pareto to solve for $\hat{a}$. Note that in the equivalent Pareto, $x_m=1$. In fact we can even take logs just and solve for the exponential scale parameter. Either way this gives that $\hat{a}$ is the mean of the $\log(z_i)$.

So any time we have an estimate of $\hat{b}$ we can immediately get the corresponding ML estimate $\hat{a}$ using step 1.

2. This leaves us with "how to estimate $b$". Note that given $c$, the discussion for estimating $x_m$ under the Pareto would lead us to set $\hat{b}=\min(x_i-c)=\min(x_i)-c=x_{(1)}-c$.

(This precludes the case where we estimate $c$ by setting $\hat{c}=\min(x_i)$ unless we then remove that observation from further participation, but I don't think you'll have that issue. I discuss at the end how to deal with that, though.)

The same thing is used when estimating location and scale for a shifted exponential; set location to the minimum observation and then use the remaining observations to estimate the scale parameter.

Note that this doesn't satisfy $z_i>1$ for the observation we just used (but does for all the others). We actually need the same trick we were just discussing, except this time it's exactly the shifted exponential issue; we eliminate the observation that gives $\hat{b}$, i.e. set $\hat{b}=x_{(1)}-c$ and lose that observation when estimating $a$.

If you estimate both $b$ and $c$ from the data you will lose the smallest two observations, the first in obtaining $c$ ($\hat{c}=x_{(1)}$) and the second obtaining $b$ ($\hat{b}=x_{(2)}-\hat{c}=x_{(2)}-x_{(1)}$), and then estimate $a$ as the mean of the logs of the remaining $z_i$.

• I will first attempt to execute the minimization problem using the partials first then if necessary look at what you just wrote. I think trying to force my problem into open source libraries like GSL, dlib and nlopt. Thanks. Nov 3, 2015 at 23:51
• What I have there should maximize the likelihood; one difficulty for the usual iterative methods is that the the ML estimate for $b$ is at a boundary (but as long as you deal with any attempt by the minimizer to cross that boundary it should still work; note that $-\log \mathcal{L}$ will go to infinity there). However, I'm not sure why you'd need to, since the two steps I outline above will be very fast and should result in what you need (you can always check by evaluating the likelihood function in a plausible region near it - something that is good practice in any case). Nov 3, 2015 at 23:58
• Excellent I will look at this. Thanks for your help. Nov 4, 2015 at 0:32