Problem für das Anpassung von Daten mit einer logarithmischen Funktion über maximale Wahrscheinlichkeitsschätzung für Fr ⇐ Python

[Anonymous]

Ich arbeite mit einem Datensatz, in dem ich die Parameter einer lognormalen Verteilung unter Verwendung der maximalen Wahrscheinlichkeitsschätzung (MLE) schätzen möchte. Die Daten stellen die beobachteten Wahrscheinlichkeiten von Ereignissen dar (wie gering, stark und kollabierten Schäden), die auf verschiedenen Ebenen einer Prädiktorvariablen auftreten, z. Peak Ground Acceleration (PGA). Insbesondere: < /p>

Bei einem Schadensniveau endet die Optimierung mit unrealistischen Parameterschätzungen und konvergiert nicht ordnungsgemäß. Die Parameterschätzungen sind extrem hoch oder unrealistisch. Probleme? Gibt es in Python Standardmethoden oder integrierte Funktionen, die dazu beitragen können, diesen Prozess zu rationalisieren? /> < /p>
Das Problem aus meinem Skript trat beim Ausführen auf. Ich habe versucht, drei Fragilitätskurven (leicht, stark und zusammengebrochen) zu tun. Die Mittelwert- und Standardabweichungswerte für geringfügige und kollabierte sind unlogisch.

Code: Select all

print(res_slight)
final_simplex: (array([[1., 1.],
[1., 1.],
[1., 1.]]), array([inf, inf, inf]))
fun: inf
message: 'Maximum number of function evaluations has been exceeded.'
nfev: 403
nit: 101
status: 1
success: False
x: array([1., 1.])

print(res_heavily)
final_simplex: (array([[0.63531661, 0.88823413],
[0.6353701 , 0.88832485],
[0.63533704, 0.88817963]]), array([203.91066709, 203.91066948, 203.91068102]))
fun: 203.91066708806784
message: 'Optimization terminated successfully.'
nfev: 60
nit: 32
status: 0
success: True
x: array([0.63531661, 0.88823413])

print(res_collapsed)
final_simplex: (array([[8.67945396e+20, 3.53694860e+01],
[7.84728875e+20, 3.53008691e+01],
[7.54836798e+20, 3.52584876e+01]]), array([183.64705108, 183.64785074, 183.6479341 ]))
fun: 183.64705108393747
message: 'Maximum number of iterations has been exceeded.'
nfev: 1823
nit: 1000
status: 2
success: False
x:  array([8.67945396e+20, 3.53694860e+01])

This is why my fragility curves like the attached figure:


My code to derive fragility curves is:

Code: Select all

from functools import partial
import numpy as np
from scipy import optimize, stats

PGA_Values = np.array([0.12, 0.14, 0.16, 0.18, 0.20, 0.22, 0.24, 0.26, 0.28, 0.30, 0.32, 0.34, 0.36, 0.38, 0.40,
0.42, 0.44, 0.46, 0.48, 0.50, 0.52, 0.54, 0.56, 0.58, 0.60, 0.62, 0.64, 0.66, 0.68, 0.70,
0.72])

im_log = np.log(PGA_Values)

#Slight Damage Data Results
slight_damage_analyses = np.array([46, 11, 18, 54, 30, 1482, 17, 73, 1082, 341, 798, 9509, 4226, 107, 7, 1,1,
3, 11, 103, 9, 27, 21, 22, 40, 22, 90, 98, 1006, 402, 1])

slight_damage_collapses = np.array([29, 5, 12, 33, 18, 1314, 10, 47, 656, 217, 611, 6751, 2734, 69, 2, 0, 0,
0, 5, 49, 5, 0, 7, 8, 13, 3, 55, 67, 342, 260, 1])

#Heavily  Damage Data Results
heavily_damage_analyses = np.array([46, 11, 18, 54, 30, 1482, 17, 73, 1082, 341, 798, 9509, 4226, 107, 7, 1,1,
3, 11, 103, 9, 27, 21, 22, 40, 22, 90, 98, 1006, 402, 1])

heavily_damage_collapses = np.array([15, 6, 3, 18, 12, 150, 5, 18, 342, 83, 160, 2244, 1107, 27, 1, 1, 1, 1,
4, 49, 3, 6, 10, 12, 17, 6, 34, 27, 597, 114, 0])

#Collapsed Data Results
collapsed_damage_analyses = np.array([46, 11, 18, 54, 30, 1482, 17, 73, 1082, 341, 798, 9509, 4226, 107, 7, 1,1,
3, 11, 103, 9, 27, 21, 22, 40, 22, 90, 98, 1006, 402, 1])

collapsed_damage_collapses = np.array([2, 0, 3, 3, 0, 18, 2, 8, 84, 41, 27, 514, 385, 11, 4, 1, 0, 2, 2, 5, 1,
21, 4, 2, 10, 13, 1, 4, 67, 28, 0])

def neg_log_likelihood_sum(params, im_l, no_a, no_c):
sigma, beta = params
theoretical_fragility_function = stats.norm(np.log(sigma), beta).cdf(im_l)
likelihood = stats.binom.pmf(no_c, no_a, theoretical_fragility_function)
log_likelihood = np.log(likelihood)
log_likelihood_sum = np.sum(log_likelihood + 1e+10)

print(FORMAT_STRING.format(sigma, beta, log_likelihood_sum))
return -log_likelihood_sum

# Slight Damage Data fitting using MLE
neg_log_likelihood_sum_partial_slight = partial(neg_log_likelihood_sum, im_l=np.log(PGA_Values), no_a=slight_damage_analyses, no_c=slight_damage_collapses)
res_slight = optimize.minimize(neg_log_likelihood_sum_partial_slight, (1, 1), method="Nelder-Mead")

# Heavily Damage Data fitting using MLE
neg_log_likelihood_sum_partial_heavily = partial(neg_log_likelihood_sum, im_l=np.log(PGA_Values), no_a=heavily_damage_analyses, no_c=heavily_damage_collapses)
res_heavily = optimize.minimize(neg_log_likelihood_sum_partial_heavily, (1, 1), method="Nelder-Mead")

# Collapsed Data fitting using MLE
neg_log_likelihood_sum_partial_collapsed = partial(neg_log_likelihood_sum, im_l=np.log(PGA_Values), no_a=collapsed_damage_analyses, no_c=collapsed_damage_collapses)
res_collapsed = optimize.minimize(neg_log_likelihood_sum_partial_collapsed, (1, 1), method="Nelder-Mead")

# Plotting Fragility Curves
import matplotlib.pyplot as plt

x = np.linspace(0, 1.0, 100)

# Slight Damage
y_slight = stats.norm(np.log(res_slight.x[0]), res_slight.x[1]).cdf(np.log(x))
plt.plot(x, y_slight, label='Slight Damage')

# Heavily Damage
y_heavily = stats.norm(np.log(res_heavily.x[0]), res_heavily.x[1]).cdf(np.log(x))
plt.plot(x, y_heavily, label='Heavily Damage')

# Collapsed
y_collapsed = stats.norm(np.log(res_collapsed.x[0]), res_collapsed.x[1]).cdf(np.log(x))
plt.plot(x, y_collapsed, label='Collapsed')

plt.scatter(PGA_Values, slight_damage_collapses / slight_damage_analyses, color='blue', alpha=0.5)
plt.scatter(PGA_Values, heavily_damage_collapses / heavily_damage_analyses, color='orange', alpha=0.5)
plt.scatter(PGA_Values, collapsed_damage_collapses / collapsed_damage_analyses, color='r', alpha=0.5)