[PYTHON] Versuchen wir es noch einmal. Wahrscheinlichste Schätzung und Anpassung des Modells (Wahrscheinlichkeitsverteilung) ① Diskrete Wahrscheinlichkeitsverteilung

Höchstwahrscheinlich Schätzmethode

Eine Möglichkeit, dies zu tun, besteht darin, eine Wahrscheinlichkeitsverteilung anzuwenden. Die Form der Wahrscheinlichkeitsverteilung ändert sich jedoch in Abhängigkeit von den Parametern, so dass es notwendig ist, die geeigneten Parameterwerte basierend auf den beobachteten Daten für die Anpassung zu bestimmen. Die wahrscheinlichste Schätzmethode ist eine Hauptmethode zur Berechnung des Parameterwerts.

Bei grober Beschreibung des Verfahrens des wahrscheinlichsten Schätzverfahrens unter der Annahme, dass die erhaltenen Daten N Daten Y sind, die aus einer Wahrscheinlichkeitsverteilung mit einem bestimmten Parameter & thgr; extrahiert wurden, beträgt die Extraktionswahrscheinlichkeit jedes Werts $ p (y_i | & thgr;). ) Die Wahrscheinlichkeitsfunktion L sei das Produkt von $ mal N Stück. $L(θ|{y_i}) = p(y_1|θ) × p(y_2|θ) × ··· × p(y_N|θ)=\prod_{i=1}^{N}p(y_N|θ)$ Berechnen Sie den Parameterwert, der die Wahrscheinlichkeitsfunktion L maximiert. Der Ablauf besteht darin, in eine Form zu konvertieren (zu logarithmieren), die einfach zu berechnen, teilweise durch die Anzahl der Parameter zu differenzieren und analytisch zu lösen ist.

Eigentlich versuchen

Abhängig von den Daten gibt es zwei Arten von Wahrscheinlichkeitsverteilungen, den diskreten Typ und den kontinuierlichen Typ. Dieses Mal führen wir als Beispiel für die diskrete Typwahrscheinlichkeitsverteilung die wahrscheinlichste Schätzung der Poisson-Verteilung und die Anpassung der Wahrscheinlichkeitsverteilung durch.

--Diskrete Wahrscheinlichkeitsverteilung: Poisson-Verteilung (diesmal)

• ~~ Kontinuierliche Wahrscheinlichkeitsverteilung: Normalverteilung ~~ (nächstes Mal)

Diskrete Wahrscheinlichkeitsverteilung: Poisson-Verteilung

Die Poisson-Verteilung wird hauptsächlich für Zähldaten verwendet. Die Wahrscheinlichkeitsfunktion, die die Verteilung darstellt, ist wie folgt. $ P (y) = \ frac {λ ^ ye ^ {-λ}} {y!} (Parameter, der den Durchschnitt von λ = y, y = 0,1,2,3 ... darstellt) $

Die Testdaten zum Schätzen und Anpassen des Modells werden im Voraus aus der Zielwahrscheinlichkeitsverteilung erstellt, damit die Antworten abgeglichen werden können.

Python

"""Zufallszahlen aus der Poisson-Verteilung"""
import numpy as np
import matplotlib.pyplot as plt

np.random.seed(seed=10)
#Parameter λ=2.4
poisson_values = np.random.poisson(lam=2.4, size=1000)
p_y, bin_edges, patches = plt.hist(poisson_values , bins=11, range=[-0.5, 10.5], normed=True)
y_k = 0.5*(bin_edges[1:] + bin_edges[:-1])
y_k #y
p_y #p(y)

Sie können die Werte und Zahlen der 1000 erhaltenen Pseudodaten sehen. Wenden wir nun die Poisson-Verteilung auf solche Daten an und schätzen ihre Parameter am wahrscheinlichsten.

Python

"""Parameterschätzung durch partielle Differenzierung der logarithmischen Wahrscheinlichkeitsfunktion"""
import sympy #Bibliothek für algebraische Berechnungen (Lösen von Differentialen und Gleichungen)
#Variablen definieren
sympy.var('λ y ')
#Haftung p(Parameter|y)Definieren
f = (λ**y)*sympy.exp(-λ)/sympy.factorial(y)
#Logistik
logf=sympy.log(f)
#Differenziere f teilweise und erweitere die Gleichung
pdff = sympy.expand(sympy.diff(logf, λ))

Die Formel pdff, die die teilweise Differenzierung der logarithmischen Wahrscheinlichkeit in Bezug auf den Parameter λ darstellt, lautet $ \ frac {y_i} {λ} -1 $. Da dies berechnet werden kann als λ, was $ \ sum_ {i = 1} ^ {1000} (\ frac {y_i} {λ} -1) = 0 $ ist, was die Summe der Daten ist.

Python

def L_sympy(f,var,values):       
    likelihood = 0 #Anfangswert der Wahrscheinlichkeit
    for i in np.arange(len(values)): #Für die Anzahl der Daten
        # model output 
        # print(values[i])
        likelihood += f.subs(var,values[i]) #Ersetzen Sie y durch einen Wert
        # print(likelihood)
    param = sympy.solve(likelihood, λ) #Löse die Gleichung für λ
    # print(param)
    return param

L_sympy(pdff,"y",poisson_values)

Parameterschätzung:[289/125]=2.312

Der Parameter λ wurde auf 2,312 geschätzt. Die Einstellung von λ dieser Pseudodaten ist 2,4, sodass Sie sie gut abschätzen können.

Wenn das statistische Modell komplex ist

Bei der Lehrbuchmethode können Sie direkt die Formel erhalten, mit der die Parameter berechnet werden können, indem Sie die teilweise differenzierte Formel (Ableitung) wie oben beschrieben lösen. Bei der allgemeinen statistischen Modellierung müssen jedoch mehrere Parameter geschätzt werden, und die Formel für die Wahrscheinlichkeitsverteilung ist komplizierter, sodass sie nicht einfach zu lösen ist.

Statistisches Modell und Wahrscheinlichkeitsverteilung

In dem allgemein verwendeten verallgemeinerten linearen Modell (GLM) wird das statistische Modell, das die Daten beschreibt, durch die Form der Wahrscheinlichkeitsverteilung und der linearen Vorhersage (und der Verknüpfungsfunktion) dargestellt.

Das obige Beispiel behandelt Y nur als Y ... Man kann sagen, dass es sich um ein statistisches Modell handelt.

Parameterschätzung, wenn das statistische Modell kompliziert ist

Übrigens ist es nicht einfach zu lösen, da als allgemeines statistisches Modell mehrere Parameter zu schätzen sind und die Formel kompliziert ist (keine lineare Formel wie das sogenannte $ \ frac {y} {λ} -1 $, sondern $ λ Nichtlinear wie ^ 2 + θ ^ 3 + ... $). Daher verwenden wir die Methode zur Lösung nichtlinearer Gleichungen.

Python

"""Funktion der diskreten Wahrscheinlichkeitsverteilung (Wahrscheinlichkeitsmassenfunktion)"""
def probability_poisson(y,λ):
    from scipy.special import factorial
    # λ:Parameter
    # y:Daten y: Tritt y-mal auf
    # return:Daten Y Wahrscheinlichkeit P.(y|Parameter) 
    return (λ**y)*np.exp(-λ)/factorial(y)

"""Log-ähnliche Wahrscheinlichkeitsfunktion: diskreter Typ"""
def L_func(param,y):       
    likelihood = 0 #Anfangswert der Wahrscheinlichkeit
    for i in np.arange(len(y)):
        # model output 
        p = probability_poisson(y[i], param)
        # likelihood
        likelihood += -np.log(p) #Haftung
    return likelihood 

"""Parameter Schätzung"""
"""
Semi-Newton-Methode (BFGS, L)-BFGS-Methode: In komplizierten Fällen kann Speicherplatz gespart werden.
"""
from scipy import optimize
x0 = [2] #Anfangswerte der Parameter
bound = [(0, None)] #Umfang der optimalen Parametersuche(min,max)

params_MLE = optimize.minimize(L_func,x0,args=(poisson_values),method='l-bfgs-b',
                       jac=None, bounds=bound, tol=None, callback=None,
                        options={'disp': None, 'maxls': 20, 'iprint': -1,
                                 'gtol': 1e-05, 'eps': 1e-08, 'maxiter': 15000,
                                 'ftol': 2.220446049250313e-09, 'maxcor': 10,
                                 'maxfun': 15000})

#Höchstwahrscheinlich Schätzung des Parameters
print('Parameter:',params_MLE.x)
#Anzahl der Parameter
k=1
#AIC (das Modell mit dem kleinsten Wert passt gut)
print('AIC：',params_MLE.fun*(2)+2*k)

Parameter:[2.31200054]
AIC： [3607.0081302]

In ähnlicher Weise wurde der Parameter λ auf 2,312 geschätzt.

Ich möchte die Parameter schnell kennen ...

Python

"""curve_Es gibt auch eine Parameterberechnung mit vorläufiger Anpassung"""
from scipy.optimize import curve_fit
parameters, cov_matrix = curve_fit(f=probability_poisson, xdata=y_k, ydata=p_y) 
print("Parameter:",parameters, "Mitverteilt:",cov_matrix)

Parameter: [2.31643586]Mitverteilt: [[0.00043928]]

Wenden Sie auf die Daten an und versuchen Sie zu veranschaulichen

Wenden wir das Modell (Poisson-Verteilung) auf die Daten mit dem geschätzten Parameter λ = 2,312 an.

Python

"""Illustration, ob es gilt"""
acc_mle = probability_poisson(y_k, params_MLE.x)
plt.hist(poisson_values , bins=11, range=[-0.5, 10.5], normed=True)
plt.plot(y_k,acc_mle)

Es fühlt sich gut an.