Example of python code for exponential distribution and maximum likelihood estimation (MLE)

Here, I would like to solve the basic problem of exponential distribution and the problem of maximum likelihood estimator related to it with python.

Question 1

python

N's homepage access interval T(Unit: time)Is an exponential distribution\\
f(t) = 2 e^{-2t}~~~~ (t>0)\\

To obey. At this time P(T\leq 3)Ask for.

First of all, various preparations.

python

import sympy as sym
from sympy.plotting import plot
sym.init_printing(use_unicode=True)
%matplotlib inline
oo = sym.oo #Infinity
(x,t) = sym.symbols('x t')

Next, let's draw a graph of the probability density function.

python

expr =2* sym.exp((-2*t))
#The graph of the obtained function is illustrated from 0 to 3.
plot(expr, (t, 0, 3))
#Also check that it is a probability density function just in case.
sym.integrate(expr, (t, 0, oo))

Solve the last given problem.

python

sym.integrate(expr, (t, 0, 3))

Problem 2

Now let's solve the maximum likelihood estimator problem.

python

N's homepage access interval T(Unit: time)Is an exponential distribution\\
f(t) = x e^{-xt}~~~~ (t>0)\\

Suppose you follow. Arbitrarily enter a set of independent observations of T and find the maximum likelihood estimator of x corresponding to that set of values.

Define the likelihood function. The variable changes to x instead of t. (Usually use λ instead of x.)

python

#First, enter the set of observed values.
A = list(map(float, input().split()))
#Next, define the likelihood function.
expr = 1
for t in A:
    expr = expr*(x * sym.exp((-x)*t))

print(expr)
#The graph of the likelihood function with x as a variable is shown from 0 to 3. The number 3 has no particular meaning here. Anything is fine.
plot(expr, (x, 0, 3))

Next, define the log-likelihood function.

python

L = sym.log(expr)
print(L)
#The graph of the obtained function is illustrated.
plot(L, (x, 0, 3))

Calculate where it differentiates to zero.

python

L_1 = L.diff(x,1)
plot(L_1, (x, 1, 3))
sym.solve(L_1, x)

Normally, you have to calculate the second derivative as follows and make sure that it is always negative. However, although this calculation is easy for humans, it does not work because Python does not transform the formula in the middle to make the calculation easier. (´ · ω · `) Shobon.

python

L_2 = L_1.diff(x,1)
plot(L_2, (x, 1, 2))
print(L_2)