I need a parameter to use in Poisson. I want to calculate this parameter using the Gauss-Newton algorithm, starting from a nonlinear system. I want to point out that i am assuming that the distribution is Poisson with parameter $\mu$

My goal is to use Poisson (which is not relevant to this question, because with this question i just want to get the only parameter of Gauss-Netwon to use in Poisson) to individually calculate the probabilities of a football club scoring 0 goals, 1 goal, 2 goals, 3 goals, 4 goals, etc. etc. So I need to find a parameter using the Gauss-Newton algorithm.

In an experimental test, i tried to use the Arithmetic Mean in Poisson, then i use the parameter 1.6 to calculate the probability that a football team individually scores goals, 1 goal, 2 goals, 3 goals, 4 goals, etc. etc. The results are decent, but I would like to get something more precise. So I thought of using another parameter, which is what i would like to find with the Gauss-Newton algorithm and the nonlinear system.

The basic, initial idea, is to start from List_Goals: [1, 2, 2, 1, 2] which i organize in the non-linear system, building it with the goals [1, 2, 2, 1, 2] and the number of times the event occurs


List_Goals: [1, 2, 2, 1, 2]
Matches played: 5 
Standard Deviation: 0.55


If a football club, in the last 5 matches,scores [1, 2, 2, 1, 2], it means that:

  • 0 goals, scored 0 times ( );
  • 1 or less goals, scored 2 times (1, 1);
  • 2 or less goals, scored 5 times (1, 2, 2, 1, 2);


If the goals scored are <= 0 events, then i will have 0/5 = 0;
If the goals scored are <= 1 events, then i will have 2/5 = 0.4;
If the goals scored are <= 2 events, then i will have 5/5 = 1;

f(0) = 0/5 = 0;
f(1) = 2/5= 0.4;
f(2) = 5/5= 1;

System: {f(0) = F(0) } therefore 0/5= 0
        {f(1) = F(1) } therefore 2/5= 0.4
        {f(2) = F(2) } therefore 5/5= 1

In this case i have three equations in one unknown, but obviously the system must be set up with a number of equations ranging from F(0) to F(max_goals_scored). This way I would have max_goals_scored + 1 equations.I should therefore start the solution starting from max_goals_scored+1 equations and increase, for example up to a maximum of 5 goals:

F(0), F(1), F(2), F(3), F(4), F(5) 


I'm having trouble using the Gauss-Newton algorithm, so i'm asking for help. I understand the theory, but i'm having trouble putting it into practice. Specifying that, what I'm looking for for the resolution and the final result is this:

a) I identify the Jacobian matrix J, which in general is composed of the partial derivatives with respect to all the parameters. In this case i only have one parameter and the matrix is actually a column vector;

b) I calculate the transpose of J multiplied by J, which in this case is a scalar (formally a number, which however in this case is a function of the unknown parameter);

c) algorithmically speaking it looks good: the inverse of J^T*J is still a scalar, so I have to divide by that number (which however will contain the mu parameter);

d) I write the iterative method according to the formula and a stopping criterion;

Naturally, these steps therefore, include functions and functional operators (the derivative). How can i use the Gauss-Network algorithm on the non-linear system and get the final result? Can you show me the various steps with my numerical data applied? (and not just with the formulas). Thank you all!!!

  • 1
    $\begingroup$ Sorry but it's not clear to me what system you are actually trying to solve? Can you write clearly the objective function you wish to minimize? $\endgroup$ Oct 30, 2023 at 13:52
  • $\begingroup$ @giorgi nguyen I have updated the question, including more information and also an image. I'm not very knowledgeable about nonlinear systems, so I apologize if I can't give you more information. I hope that what I wrote in UPDATE is useful for you to understand what I would like to achieve $\endgroup$ Oct 30, 2023 at 14:56
  • $\begingroup$ @giorginguyen I updated the question again. I hope you can help me $\endgroup$ Oct 31, 2023 at 5:39

1 Answer 1


This is not an answer, but is too long for a comment.

It is not entirely clear to me what you are asking, you seem to have all the information necessary. However the notation could be improved.

I am using $F_k(\mu) = e^{-\mu} \sum_{j=0}^k { \mu^j \over j!}$ and am assuming that $f_e(k)$ represents the corresponding experimental datum. It is straightforward to compute $F_k'(\mu)$.

I am guessing that you have a system of equations $r:\mathbb{R} \to \mathbb{R}^3$, $r_k(\mu) = F_k(\mu)-f_e(k)$, then you would like to solve $r(\mu) \approx 0$.

To apply at iteration of Gauss Newton at some iteration $e_n$ you linearise $r$ as in $r(\mu_{n+1}) \approx r(\mu_n) + Dr(\mu_n) (\mu_{n+1}-\mu_n)$ and (attempt to) solve the system $r(\mu_n) + Df(\mu_n) (\mu_{n+1}-\mu_n) = 0$ to compute $e_{n+1}$.

Using ordinary least squares, we have the update rule $\mu_{n+1} = \mu_n - (Dr(\mu_n)^T Dr(\mu_n))^{-1} Dr(\mu_n)^T r(\mu_n)$.

In this case (assuming that I understand correctly) you have, $Dr(\mu) = \begin{bmatrix} F_0'(\mu) \\ F_1'(\mu) \\ F_2'(\mu) \end{bmatrix}$.

A readable explanation of Gauss Newton can be found here, some standard improvements (Armijo step size, LevenbergMarquardt) can be found here.

Another approach:

The following plots the sum of errors squared vs. the Poisson $\mu$ parameter for the data [1, 2, 2, 1, 2]. There could easily be errors in my implementation.

#!/usr/bin/env python3

import matplotlib.pyplot as plt
import numpy as np
import math

goals = [1, 2, 2, 1, 2]
def f(k):
    """Observed data."""
    return len(list(filter(lambda j: j <= k, goals)))/len(goals)

def poisson_cdf(mu, k):
    """Poisson CDF."""
    return sum(mu**j/math.factorial(j) for j in range(0, k+1))*math.exp(-mu)

def lsq_error(mu, k_vals):
    """Sum of squares of errors."""
    return sum((poisson_cdf(mu, k)-f(k))**2 for k in k_vals)

# Find the error for values of k in k_vals.
k_vals = [0, 1, 2]
# Plot x range.
mu_vals = np.arange(0.0, 5.0, 0.1)
lsq_error_vals = np.array(list(lsq_error(mu, k_vals) for mu in mu_vals))

fig, ax = plt.subplots()
ax.plot(mu_vals, lsq_error_vals)

ax.set(xlabel='mu', ylabel='lsq error',
       title='sum of errors squared vs. mu')


This gives the following plot from which we can visually see that there is a $\min$ around $1.6$. enter image description here

  • $\begingroup$ Let us continue this discussion in chat. $\endgroup$ Nov 2, 2023 at 16:40
  • 3
    $\begingroup$ @Horiatiki copper.hat attempted to help you. It seems that their best efforts have not been to your liking, which is fine. They have given you as much help as they are willing and/or able to do. At this point, you need to stop pestering them, as I think they have made it clear that they no longer wish to engage with this question. $\endgroup$
    – Xander Henderson
    Nov 2, 2023 at 19:24

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .