Open In Colab

Berdasarkan isu #126: anfrek: Log Pearson III

Referensi Isu:

  • Soetopo, W., Montarcih, L., Press, U. B., & Media, U. (2017). Rekayasa Statistika untuk Teknik Pengairan. Universitas Brawijaya Press. https://books.google.co.id/books?id=TzVTDwAAQBAJ
  • Limantara, L. (2018). Rekayasa Hidrologi.
  • Soewarno. (1995). hidrologi: Aplikasi Metode Statistik untuk Analisa Data. NOVA.

Deskripsi Isu:

  • Mencari nilai ekstrim dengan kala ulang tertentu menggunakan distribusi Log Pearson III. Penerapan ini bisa digunakan untuk hujan rancangan atau debit banjir rancangan.

Diskusi Isu:

  • #156 - Bagaimana menghitung periode ulang distribusi (analisis frekuensi) tanpa melihat tabel?

Strategi:

  • Luaran dari fungsi merupakan tabel atau array dengan kala ulang yang dapat dijadikan sebagai input, jika tidak menggunakan kala ulang yang umum (2, 5, ..., 100)
  • Buat hasil perhitungan berdasarkan manual (tabel) dan menggunakan fungsi yang tersedia (scipy).

PERSIAPAN DAN DATASET

In [1]:
import numpy as np
import pandas as pd
from scipy import stats
In [2]:
# contoh data diambil dari buku
# Rekayasa Hidrologi hal. 109

_HUJAN = np.array([85, 92, 115, 116, 122, 52, 69, 95, 96, 105])
_TAHUN = np.arange(1998, 2008) # 1998-2007

data = pd.DataFrame(
    data=np.stack([_TAHUN, _HUJAN], axis=1),
    columns=['tahun', 'hujan']
)
data.tahun = pd.to_datetime(data.tahun, format='%Y')
data.set_index('tahun', inplace=True)
data
Out[2]:
hujan
tahun
1998-01-01 85
1999-01-01 92
2000-01-01 115
2001-01-01 116
2002-01-01 122
2003-01-01 52
2004-01-01 69
2005-01-01 95
2006-01-01 96
2007-01-01 105

TABEL

Terdapat 2 tabel untuk modul hk126 yaitu:

  • t_pearson3_sw: Tabel Nilai $K$ Distribusi Pearson Tipe III dari buku hidrologi: Aplikasi Metode Statistik untuk Analisa Data oleh Soewarno.
  • t_pearson3_st: Tabel Distribusi Pearson Type III (nilai $K$) dari buku Rekayasa Statistika untuk Teknik Pengairan.
  • t_pearson3_lm: Tabel Distribusi Log Pearson Tipe III Nilai G untuk Cs Positif dan Negatif

Dalam modul hk126 nilai $K$ akan dibangkitkan menggunakan fungsi scipy.stats.pearson3.ppf secara default. Mohon diperhatikan jika ingin menggunakan nilai $K$ yang berasal dari sumber lain.

In [3]:
# tabel dari Soewarno, hal. 219
# lampiran tabel III-3, Nilai k Distribusi Pearson tipe III
# dan Log Pearson ti

# kode: SW

_DATA_SW = np.array([
    [-0.360, 0.420, 1.180, 2.278, 3.152, 4.051, 4.970, 7.250],
    [-0.360, 0.518, 1.250, 2.262, 3.048, 3.845, 4.652, 6.600],
    [-0.330, 0.574, 1.284, 2.240, 2.970, 3.705, 4.444, 6.200],
    [-0.307, 0.609, 1.302, 2.219, 2.912, 3.605, 4.298, 5.910],
    [-0.282, 0.643, 1.318, 2.193, 2.848, 3.499, 4.147, 5.660],
    [-0.254, 0.675, 1.329, 2.163, 2.780, 3.388, 3.990, 5.390],
    [-0.225, 0.705, 1.337, 2.128, 2.706, 3.271, 3.828, 5.110],
    [-0.195, 0.732, 1.340, 2.087, 2.626, 3.149, 3.661, 4.820],
    [-0.164, 0.758, 1.340, 2.043, 2.542, 3.022, 3.489, 4.540],
    [-0.148, 0.769, 1.339, 2.018, 2.498, 2.957, 3.401, 4.395],
    [-0.132, 0.780, 1.336, 1.998, 2.453, 2.891, 3.312, 4.250],
    [-0.116, 0.790, 1.333, 1.967, 2.407, 2.824, 3.223, 4.105],
    [0.099, 0.800, 1.328, 1.939, 2.359, 2.755, 3.132, 3.960],
    [-0.083, 0.808, 1.323, 1.910, 2.311, 2.686, 3.041, 3.815],
    [-0.066, 0.816, 1.317, 1.880, 2.261, 2.615, 2.949, 3.670],
    [-0.050, 0.824, 1.309, 1.849, 2.211, 2.544, 2.856, 3.525],
    [-0.033, 0.830, 1.301, 1.818, 2.159, 2.472, 2.763, 3.380],
    [-0.017, 0.836, 1.292, 1.785, 2.107, 2.400, 2.670, 3.235],
    [0.000, 0.842, 1.282, 1.751, 2.054, 2.326, 2.576, 3.090],
    [0.017, 0.836, 1.270, 1.761, 2.000, 2.252, 2.482, 3.950],
    [0.033, 0.850, 1.258, 1.680, 1.945, 2.178, 2.388, 2.810],
    [0.050, 0.853, 1.245, 1.643, 1.890, 2.104, 2.294, 2.675],
    [0.066, 0.855, 1.231, 1.606, 1.834, 2.029, 2.201, 2.540],
    [0.083, 0.856, 1.216, 1.567, 1.777, 1.955, 2.108, 2.400],
    [0.099, 0.857, 1.200, 1.528, 1.720, 1.880, 2.016, 2.275],
    [0.116, 0.857, 1.183, 1.488, 1.663, 1.806, 1.926, 2.150],
    [0.132, 0.856, 1.166, 1.448, 1.606, 1.733, 1.837, 2.035],
    [0.148, 0.854, 1.147, 1.407, 1.549, 1.660, 1.749, 1.910],
    [0.164, 0.852, 1.128, 1.366, 1.492, 1.588, 1.664, 1.800],
    [0.195, 0.844, 1.086, 1.282, 1.379, 1.449, 1.501, 1.625],
    [0.225, 0.832, 1.041, 1.198, 1.270, 1.318, 1.351, 1.465],
    [0.254, 0.817, 0.994, 1.116, 1.166, 1.197, 1.216, 1.280],
    [0.282, 0.799, 0.945, 1.035, 1.069, 1.087, 1.097, 1.130],
    [0.307, 0.777, 0.895, 0.959, 0.980, 0.990, 1.995, 1.000],
    [0.330, 0.752, 0.844, 0.888, 0.900, 0.905, 0.907, 0.910],
    [0.360, 0.711, 0.771, 0.793, 0.798, 0.799, 0.800, 0.802],
    [0.396, 0.636, 0.660, 0.666, 0.666, 0.667, 0.667, 0.668]]
)

_INDEX_SW = [
    3, 2.5, 2.2, 2, 1.8, 1.6, 1.4, 1.2, 1,
    0.9, 0.8, 0.7, 0.6, 0.5, 0.4, 0.3, 0.2, 0.1, 0,
    -0.1, -0.2, -0.3, -0.4, -0.5, -0.6, -0.7, -0.8, -0.9, -1,
    -1.2, -1.4, -1.6, -1.8, -2. , -2.2, -2.5, -3.
]

_COL_SW = [0.5, 0.2, 0.1, 0.04, 0.02, 0.01, 0.005, 0.001]

t_pearson3_sw = pd.DataFrame(data=_DATA_SW, index=_INDEX_SW, columns=_COL_SW)
t_pearson3_sw
Out[3]:
0.500 0.200 0.100 0.040 0.020 0.010 0.005 0.001
3.0 -0.360 0.420 1.180 2.278 3.152 4.051 4.970 7.250
2.5 -0.360 0.518 1.250 2.262 3.048 3.845 4.652 6.600
2.2 -0.330 0.574 1.284 2.240 2.970 3.705 4.444 6.200
2.0 -0.307 0.609 1.302 2.219 2.912 3.605 4.298 5.910
1.8 -0.282 0.643 1.318 2.193 2.848 3.499 4.147 5.660
1.6 -0.254 0.675 1.329 2.163 2.780 3.388 3.990 5.390
1.4 -0.225 0.705 1.337 2.128 2.706 3.271 3.828 5.110
1.2 -0.195 0.732 1.340 2.087 2.626 3.149 3.661 4.820
1.0 -0.164 0.758 1.340 2.043 2.542 3.022 3.489 4.540
0.9 -0.148 0.769 1.339 2.018 2.498 2.957 3.401 4.395
0.8 -0.132 0.780 1.336 1.998 2.453 2.891 3.312 4.250
0.7 -0.116 0.790 1.333 1.967 2.407 2.824 3.223 4.105
0.6 0.099 0.800 1.328 1.939 2.359 2.755 3.132 3.960
0.5 -0.083 0.808 1.323 1.910 2.311 2.686 3.041 3.815
0.4 -0.066 0.816 1.317 1.880 2.261 2.615 2.949 3.670
0.3 -0.050 0.824 1.309 1.849 2.211 2.544 2.856 3.525
0.2 -0.033 0.830 1.301 1.818 2.159 2.472 2.763 3.380
0.1 -0.017 0.836 1.292 1.785 2.107 2.400 2.670 3.235
0.0 0.000 0.842 1.282 1.751 2.054 2.326 2.576 3.090
-0.1 0.017 0.836 1.270 1.761 2.000 2.252 2.482 3.950
-0.2 0.033 0.850 1.258 1.680 1.945 2.178 2.388 2.810
-0.3 0.050 0.853 1.245 1.643 1.890 2.104 2.294 2.675
-0.4 0.066 0.855 1.231 1.606 1.834 2.029 2.201 2.540
-0.5 0.083 0.856 1.216 1.567 1.777 1.955 2.108 2.400
-0.6 0.099 0.857 1.200 1.528 1.720 1.880 2.016 2.275
-0.7 0.116 0.857 1.183 1.488 1.663 1.806 1.926 2.150
-0.8 0.132 0.856 1.166 1.448 1.606 1.733 1.837 2.035
-0.9 0.148 0.854 1.147 1.407 1.549 1.660 1.749 1.910
-1.0 0.164 0.852 1.128 1.366 1.492 1.588 1.664 1.800
-1.2 0.195 0.844 1.086 1.282 1.379 1.449 1.501 1.625
-1.4 0.225 0.832 1.041 1.198 1.270 1.318 1.351 1.465
-1.6 0.254 0.817 0.994 1.116 1.166 1.197 1.216 1.280
-1.8 0.282 0.799 0.945 1.035 1.069 1.087 1.097 1.130
-2.0 0.307 0.777 0.895 0.959 0.980 0.990 1.995 1.000
-2.2 0.330 0.752 0.844 0.888 0.900 0.905 0.907 0.910
-2.5 0.360 0.711 0.771 0.793 0.798 0.799 0.800 0.802
-3.0 0.396 0.636 0.660 0.666 0.666 0.667 0.667 0.668
In [4]:
# Tabel dari Soetopo hal. 105
# Tabel Distribusi Pearson Type III (nilai K)

# KODE: ST

_DATA_ST = [
    [-4.051, -2.003, -1.180, -0.420, 0.396, 0.636, 0.660, 0.666, 0.666, 0.667, 0.667, 0.667],
    [-4.013, -2.007, -1.195, -0.440, 0.390, 0.651, 0.681, 0.688, 0.689, 0.690, 0.690, 0.690],
    [-3.973, -2.010, -1.210, -0.460, 0.384, 0.666, 0.702, 0.712, 0.714, 0.714, 0.714, 0.714],
    [-3.932, -2.012, -1.224, -0.479, 0.376, 0.681, 0.724, 0.738, 0.740, 0.740, 0.741, 0.741],
    [-3.889, -2.013, -1.238, -0.499, 0.369, 0.696, 0.747, 0.765, 0.768, 0.769, 0.769, 0.769],
    [-3.845, -2.012, -1.250, -0.518, 0.360, 0.711, 0.771, 0.793, 0.798, 0.799, 0.800, 0.800],
    [-3.800, -2.011, -1.262, -0.537, 0.351, 0.725, 0.795, 0.823, 0.830, 0.832, 0.833, 0.833],
    [-3.753, -2.009, -1.274, -0.555, 0.341, 0.739, 0.819, 0.855, 0.864, 0.867, 0.869, 0.869],
    [-3.705, -2.006, -1.284, -0.574, 0.330, 0.752, 0.844, 0.888, 0.900, 0.905, 0.907, 0.909],
    [-3.656, -2.001, -1.294, -0.592, 0.319, 0.765, 0.869, 0.923, 0.939, 0.946, 0.949, 0.952],
    [-3.605, -1.996, -1.303, -0.609, 0.307, 0.777, 0.895, 0.959, 0.980, 0.990, 0.995, 0.999],
    [-3.553, -1.989, -1.311, -0.627, 0.294, 0.788, 0.920, 0.997, 1.023, 1.037, 1.044, 1.051],
    [-3.499, -1.981, -1.318, -0.643, 0.282, 0.799, 0.945, 1.035, 1.069, 1.087, 1.097, 1.107],
    [-3.444, -1.972, -1.324, -0.660, 0.268, 0.808, 0.970, 1.075, 1.116, 1.140, 1.155, 1.170],
    [-3.388, -1.962, -1.329, -0.675, 0.254, 0.817, 0.994, 1.116, 1.166, 1.197, 1.216, 1.238],
    [-3.330, -1.951, -1.333, -0.691, 0.240, 0.825, 1.018, 1.157, 1.217, 1.256, 1.282, 1.313],
    [-3.271, -1.938, -1.337, -0.705, 0.225, 0.832, 1.041, 1.198, 1.270, 1.318, 1.351, 1.394],
    [-3.211, -1.925, -1.339, -0.719, 0.210, 0.838, 1.064, 1.240, 1.324, 1.383, 1.424, 1.482],
    [-3.149, -1.910, -1.340, -0.733, 0.195, 0.844, 1.086, 1.282, 1.379, 1.449, 1.501, 1.577],
    [-3.087, -1.894, -1.341, -0.745, 0.180, 0.848, 1.107, 1.324, 1.435, 1.518, 1.581, 1.678],
    [-3.023, -1.877, -1.340, -0.758, 0.164, 0.852, 1.128, 1.366, 1.492, 1.588, 1.664, 1.786],
    [-2.957, -1.859, -1.339, -0.769, 0.148, 0.854, 1.147, 1.407, 1.549, 1.660, 1.749, 1.899],
    [-2.891, -1.839, -1.336, -0.780, 0.132, 0.856, 1.166, 1.448, 1.606, 1.733, 1.837, 2.017],
    [-2.824, -1.819, -1.333, -0.790, 0.116, 0.857, 1.183, 1.489, 1.663, 1.806, 1.926, 2.141],
    [-2.755, -1.797, -1.329, -0.800, 0.099, 0.857, 1.200, 1.528, 1.720, 1.880, 2.016, 2.268],
    [-2.686, -1.774, -1.323, -0.808, 0.083, 0.857, 1.216, 1.567, 1.777, 1.955, 2.108, 2.399],
    [-2.615, -1.750, -1.317, -0.816, 0.067, 0.855, 1.231, 1.606, 1.834, 2.029, 2.201, 2.533],
    [-2.544, -1.726, -1.309, -0.824, 0.050, 0.853, 1.245, 1.643, 1.890, 2.104, 2.294, 2.669],
    [-2.472, -1.700, -1.301, -0.830, 0.033, 0.850, 1.258, 1.680, 1.945, 2.178, 2.388, 2.808],
    [-2.400, -1.673, -1.292, -0.836, 0.017, 0.846, 1.270, 1.716, 2.000, 2.253, 2.482, 2.948],
    [-2.326, -1.645, -1.282, -0.842, 0.000, 0.842, 1.282, 1.751, 2.054, 2.326, 2.576, 3.090],
    [-2.253, -1.616, -1.270, -0.846, -0.017, 0.836, 1.292, 1.785, 2.107, 2.400, 2.670, 3.233],
    [-2.178, -1.586, -1.258, -0.850, -0.033, 0.830, 1.301, 1.818, 2.159, 2.472, 2.763, 3.377],
    [-2.104, -1.555, -1.245, -0.853, -0.050, 0.824, 1.309, 1.849, 2.211, 2.544, 2.856, 3.521],
    [-2.029, -1.524, -1.231, -0.855, -0.067, 0.816, 1.317, 1.880, 2.261, 2.615, 2.949, 3.666],
    [-1.955, -1.491, -1.216, -0.857, -0.083, 0.808, 1.323, 1.910, 2.311, 2.686, 3.041, 3.811],
    [-1.880, -1.458, -1.200, -0.857, -0.099, 0.800, 1.329, 1.939, 2.359, 2.755, 3.132, 3.956],
    [-1.806, -1.423, -1.183, -0.857, -0.116, 0.790, 1.333, 1.967, 2.407, 2.824, 3.223, 4.100],
    [-1.733, -1.389, -1.166, -0.856, -0.132, 0.780, 1.336, 1.993, 2.453, 2.891, 3.312, 4.244],
    [-1.660, -1.353, -1.147, -0.854, -0.148, 0.769, 1.339, 2.018, 2.498, 2.957, 3.401, 4.388],
    [-1.588, -1.317, -1.128, -0.852, -0.164, 0.758, 1.340, 2.043, 2.542, 3.023, 3.489, 4.531],
    [-1.518, -1.280, -1.107, -0.848, -0.180, 0.745, 1.341, 2.066, 2.585, 3.087, 3.575, 4.673],
    [-1.449, -1.243, -1.086, -0.844, -0.195, 0.733, 1.340, 2.088, 2.626, 3.149, 3.661, 4.815],
    [-1.383, -1.206, -1.064, -0.838, -0.210, 0.719, 1.339, 2.108, 2.667, 3.211, 3.745, 4.955],
    [-1.318, -1.168, -1.041, -0.832, -0.225, 0.705, 1.337, 2.128, 2.706, 3.271, 3.828, 5.095],
    [-1.256, -1.131, -1.018, -0.825, -0.240, 0.691, 1.333, 2.146, 2.743, 3.330, 3.910, 5.234],
    [-1.197, -1.093, -0.994, -0.817, -0.254, 0.675, 1.329, 2.163, 2.780, 3.388, 3.990, 5.371],
    [-1.140, -1.056, -0.970, -0.808, -0.268, 0.660, 1.324, 2.179, 2.815, 3.444, 4.069, 5.507],
    [-1.087, -1.020, -0.945, -0.799, -0.282, 0.643, 1.318, 2.193, 2.848, 3.499, 4.147, 5.642],
    [-1.037, -0.984, -0.920, -0.788, -0.294, 0.627, 1.311, 2.207, 2.881, 3.553, 4.223, 5.775],
    [-0.990, -0.949, -0.895, -0.777, -0.307, 0.609, 1.303, 2.219, 2.912, 3.605, 4.298, 5.908],
    [-0.946, -0.915, -0.869, -0.765, -0.319, 0.592, 1.294, 2.230, 2.942, 3.656, 4.372, 6.039],
    [-0.905, -0.882, -0.844, -0.752, -0.330, 0.574, 1.284, 2.240, 2.970, 3.705, 4.444, 6.168],
    [-0.867, -0.850, -0.819, -0.739, -0.341, 0.555, 1.274, 2.248, 2.997, 3.753, 4.515, 6.296],
    [-0.832, -0.819, -0.795, -0.725, -0.351, 0.537, 1.262, 2.256, 3.023, 3.800, 4.584, 6.423],
    [-0.799, -0.790, -0.771, -0.711, -0.360, 0.518, 1.250, 2.262, 3.048, 3.845, 4.652, 6.548],
    [-0.769, -0.762, -0.747, -0.696, -0.369, 0.499, 1.238, 2.267, 3.071, 3.889, 4.718, 6.672],
    [-0.740, -0.736, -0.724, -0.681, -0.376, 0.479, 1.224, 2.272, 3.093, 3.932, 4.783, 6.794],
    [-0.714, -0.711, -0.702, -0.666, -0.384, 0.460, 1.210, 2.275, 3.114, 3.973, 4.847, 6.915],
    [-0.690, -0.688, -0.681, -0.651, -0.390, 0.440, 1.195, 2.277, 3.134, 4.013, 4.909, 7.034],
    [-0.667, -0.665, -0.660, -0.636, -0.396, 0.420, 1.180, 2.278, 3.152, 4.051, 4.970, 7.152]
]

_INDEX_ST = np.arange(-30, 31, 1)/10

_COL_ST = np.array([99, 95, 90, 80, 50, 20, 10, 4, 2, 1, 0.5, 0.1])/100

t_pearson3_st = pd.DataFrame(data=_DATA_ST, index=_INDEX_ST, columns=_COL_ST)
t_pearson3_st
Out[4]:
0.990 0.950 0.900 0.800 0.500 0.200 0.100 0.040 0.020 0.010 0.005 0.001
-3.0 -4.051 -2.003 -1.180 -0.420 0.396 0.636 0.660 0.666 0.666 0.667 0.667 0.667
-2.9 -4.013 -2.007 -1.195 -0.440 0.390 0.651 0.681 0.688 0.689 0.690 0.690 0.690
-2.8 -3.973 -2.010 -1.210 -0.460 0.384 0.666 0.702 0.712 0.714 0.714 0.714 0.714
-2.7 -3.932 -2.012 -1.224 -0.479 0.376 0.681 0.724 0.738 0.740 0.740 0.741 0.741
-2.6 -3.889 -2.013 -1.238 -0.499 0.369 0.696 0.747 0.765 0.768 0.769 0.769 0.769
... ... ... ... ... ... ... ... ... ... ... ... ...
2.6 -0.769 -0.762 -0.747 -0.696 -0.369 0.499 1.238 2.267 3.071 3.889 4.718 6.672
2.7 -0.740 -0.736 -0.724 -0.681 -0.376 0.479 1.224 2.272 3.093 3.932 4.783 6.794
2.8 -0.714 -0.711 -0.702 -0.666 -0.384 0.460 1.210 2.275 3.114 3.973 4.847 6.915
2.9 -0.690 -0.688 -0.681 -0.651 -0.390 0.440 1.195 2.277 3.134 4.013 4.909 7.034
3.0 -0.667 -0.665 -0.660 -0.636 -0.396 0.420 1.180 2.278 3.152 4.051 4.970 7.152

61 rows × 12 columns

In [5]:
# Dari buku Limantara hal. 107-109
# Tabel Distribsi log Pearson Tipe III Nilai G
# Untuk Cs Positif & Negatif

# KODE: LM

_DATA_LM = [
    [-0.667, -0.665, -0.660, -0.636, -0.396, 0.420, 1.180, 2.278, 3.152, 4.061, 4.970],
    [-0.690, -0.688, -0.681, -0.651, -0.390, 0.440, 1.196, 2.277, 3.134, 4.013, 4.909],
    [-0.714, -0.711, -0.702, -0.666, -0.384, 0.460, 1.210, 2.275, 3.114, 3.973, 4.847],
    [-0.740, -0.736, -0.724, -0.681, -0.376, 0.479, 1.224, 2.272, 3.097, 3.932, 4.783],
    [-0.769, -0.762, -0.747, -0.695, -0.368, 0.499, 1.238, 2.267, 3.071, 3.889, 4.718],
    [-0.799, -0.790, -0.771, -0.711, -0.360, 0.518, 1.250, 2.262, 3.048, 3.845, 4.652],
    [-0.832, -0.819, -0.795, -0.725, -0.351, 0.537, 1.262, 2.256, 3.029, 3.800, 4.584],
    [-0.867, -0.850, -0.819, -0.739, -0.341, 0.555, 1.274, 2.248, 2.997, 3.753, 4.515],
    [-0.905, -0.882, -0.844, -0.752, -0.330, 0.574, 1.284, 2.240, 2.970, 3.705, 4.454],
    [-0.946, -0.914, -0.869, -0.765, -0.319, 0.592, 1.294, 2.230, 2.942, 3.656, 4.372],
    [-0.990, -0.949, -0.896, -0.777, -0.307, 0.609, 1.302, 2.219, 2.912, 3.605, 4.298],
    [-1.037, -0.984, -0.920, -0.788, -0.294, 0.627, 1.310, 2.207, 2.881, 3.553, 4.223],
    [-1.087, -1.020, -0.945, -0.799, -0.282, 0.643, 1.318, 2.193, 2.848, 3.499, 4.147],
    [-1.140, -1.056, -0.970, -0.808, -0.268, 0.660, 1.324, 2.179, 2.815, 3.444, 4.069],
    [-1.197, -1.093, -0.994, -0.817, -0.254, 0.675, 1.329, 2.163, 2.780, 3.388, 3.990],
    [-1.256, -1.131, -1.018, -0.825, -0.240, 0.690, 1.333, 2.146, 2.745, 3.330, 3.910],
    [-1.318, -1.163, -1.041, -0.832, -0.225, 0.705, 1.337, 2.128, 2.706, 3.271, 3.828],
    [-1.388, -1.206, -1.064, -0.838, -0.210, 0.719, 1.339, 2.108, 2.666, 3.211, 3.745],
    [-1.449, -1.243, -1.086, -0.844, -0.195, 0.732, 1.340, 2.087, 2.626, 3.149, 3.661],
    [-1.518, -1.280, -1.107, -0.848, -0.180, 0.745, 1.341, 2.066, 2.585, 3.087, 3.575],
    [-1.588, -1.317, -1.128, -0.852, -0.164, 0.758, 1.340, 2.043, 2.542, 3.022, 3.489],
    [-1.660, -1.353, -1.147, -0.854, -0.148, 0.769, 1.339, 2.018, 2.498, 2.967, 3.401],
    [-1.733, -1.388, -1.166, -0.856, -0.132, 0.780, 1.336, 1.993, 2.453, 2.891, 3.312],
    [-1.806, -1.423, -1.183, -0.857, -0.116, 0.790, 1.333, 1.967, 2.407, 2.824, 3.223],
    [-1.880, -1.458, -1.200, -0.857, -0.099, 0.800, 1.328, 1.939, 2.359, 2.755, 3.123],
    [-1.965, -1.491, -1.216, -0.856, -0.083, 0.808, 1.323, 1.910, 2.311, 2.686, 3.041],
    [-2.029, -1.524, -1.231, -0.855, -0.066, 0.816, 1.317, 1.880, 2.261, 2.615, 2.949],
    [-2.104, -1.555, -1.245, -0.853, -0.050, 0.824, 1.309, 1.849, 2.211, 2.544, 2.856],
    [-2.175, -1.586, -1.258, -0.850, -0.033, 0.830, 1.301, 1.818, 2.159, 2.472, 2.763],
    [-2.225, -1.616, -1.270, -0.846, -0.017, 0.836, 1.292, 1.785, 2.107, 2.400, 2.670],
    [-2.326, -1.645, -1.282, -0.842, 0.000, 0.842, 1.282, 1.751, 2.064, 2.064, 2.576],
    [-2.400, -1.673, -1.292, -0.836, 0.017, 0.846, 1.270, 1.716, 2.000, 2.252, 2.482],
    [-2.472, -1.700, -1.301, -0.830, 0.033, 0.850, 1.258, 1.680, 1.945, 2.178, 2.388],
    [-2.544, -1.762, -1.309, -0.824, 0.050, 0.853, 1.245, 0.163, 1.890, 2.104, 2.294],
    [-2.615, -1.750, -1.317, -0.816, 0.066, 0.855, 1.231, 1.606, 1.834, 2.029, 2.201],
    [-2.686, -1.774, -1.323, -0.808, 0.083, 0.856, 1.216, 1.567, 1.777, 1.955, 2.108],
    [-2.755, -1.797, -1.328, -0.800, 0.099, 0.857, 1.200, 1.528, 1.720, 1.880, 2.016],
    [-2.824, -1.819, -1.333, -0.790, 0.116, 0.857, 1.183, 1.488, 1.633, 1.800, 1.936],
    [-2.891, -1.839, -1.336, -0.780, 0.132, 0.856, 1.166, 1.484, 1.608, 1.733, 1.837],
    [-2.957, -1.858, -1.339, -0.769, 0.148, 0.854, 1.147, 1.407, 1.549, 1.660, 1.749],
    [-3.022, -1.877, -1.340, -0.758, 0.164, 0.852, 1.108, 1.366, 1.492, 1.588, 1.664],
    [-3.087, -1.894, -1.341, -0.745, 0.180, 0.848, 1.107, 1.324, 1.435, 1.518, 1.581],
    [-3.149, -1.910, -1.340, -0.732, 0.195, 0.844, 1.086, 1.282, 1.379, 1.449, 1.501],
    [-3.211, -1.925, -1.339, -0.719, 0.210, 0.838, 1.064, 1.240, 1.324, 1.383, 1.424],
    [-3.271, -1.938, -1.337, -0.705, 0.225, 0.832, 1.041, 1.196, 1.270, 1.316, 1.351],
    [-3.330, -1.961, -1.333, -0.690, 0.240, 0.825, 1.018, 1.157, 1.217, 1.256, 1.282],
    [-3.388, -1.962, -1.329, -0.675, 0.254, 0.817, 0.994, 1.116, 1.168, 1.197, 1.216],
    [-3.444, -1.972, -1.324, -0.660, 0.268, 0.808, 0.970, 1.075, 1.116, 1.140, 1.155],
    [-3.499, -1.981, -1.318, -0.643, 0.282, 0.799, 0.945, 1.035, 1.069, 1.087, 1.097],
    [-3.533, -1.989, -1.310, -0.627, 0.294, 0.788, 0.920, 0.996, 1.023, 1.037, 1.044],
    [-3.605, -1.996, -1.302, -0.609, 0.307, 0.777, 0.895, 0.969, 0.980, 0.990, 0.995],
    [-3.656, -2.001, -1.294, -0.592, 0.319, 0.765, 0.869, 0.923, 0.939, 0.346, 0.949],
    [-3.705, -2.006, -1.284, -0.574, 0.330, 0.732, 0.849, 0.888, 0.900, 0.905, 0.907],
    [-3.753, -2.009, -1.274, -0.555, 0.341, 0.739, 0.819, 0.855, 0.864, 0.867, 0.869],
    [-3.800, -2.011, -1.262, -0.537, 0.351, 0.725, 0.795, 0.823, 0.830, 0.832, 0.833],
    [-3.845, -2.012, -1.250, -0.518, 0.360, 0.711, 0.771, 0.793, 0.796, 0.799, 0.800],
    [-3.889, -2.013, -1.238, -0.499, 0.368, 0.696, 0.747, 0.764, 0.767, 0.769, 0.769],
    [-3.932, -2.011, -1.224, -0.479, 0.376, 0.681, 0.724, 0.738, 0.740, 0.740, 0.741],
    [-3.973, -2.010, -1.210, -0.460, 0.384, 0.666, 0.702, 0.712, 0.714, 0.734, 0.714],
    [-4.013, -2.007, -1.195, -0.440, 0.330, 0.651, 0.681, 0.683, 0.689, 0.690, 0.690],
    [-4.051, -2.003, -1.180, -0.420, 0.390, 0.636, 0.660, 0.666, 0.666, 0.667, 0.667]
]

_INDEX_LM = np.arange(30, -31, -1) / 10

_COL_LM =  np.array([99, 95, 90, 80, 50, 20, 10, 4, 2, 1, 0.5])/100

t_pearson3_lm = pd.DataFrame(data=_DATA_LM, index=_INDEX_LM, columns=_COL_LM)
t_pearson3_lm
Out[5]:
0.990 0.950 0.900 0.800 0.500 0.200 0.100 0.040 0.020 0.010 0.005
3.0 -0.667 -0.665 -0.660 -0.636 -0.396 0.420 1.180 2.278 3.152 4.061 4.970
2.9 -0.690 -0.688 -0.681 -0.651 -0.390 0.440 1.196 2.277 3.134 4.013 4.909
2.8 -0.714 -0.711 -0.702 -0.666 -0.384 0.460 1.210 2.275 3.114 3.973 4.847
2.7 -0.740 -0.736 -0.724 -0.681 -0.376 0.479 1.224 2.272 3.097 3.932 4.783
2.6 -0.769 -0.762 -0.747 -0.695 -0.368 0.499 1.238 2.267 3.071 3.889 4.718
... ... ... ... ... ... ... ... ... ... ... ...
-2.6 -3.889 -2.013 -1.238 -0.499 0.368 0.696 0.747 0.764 0.767 0.769 0.769
-2.7 -3.932 -2.011 -1.224 -0.479 0.376 0.681 0.724 0.738 0.740 0.740 0.741
-2.8 -3.973 -2.010 -1.210 -0.460 0.384 0.666 0.702 0.712 0.714 0.734 0.714
-2.9 -4.013 -2.007 -1.195 -0.440 0.330 0.651 0.681 0.683 0.689 0.690 0.690
-3.0 -4.051 -2.003 -1.180 -0.420 0.390 0.636 0.660 0.666 0.666 0.667 0.667

61 rows × 11 columns

KODE

In [6]:
# KODE FUNGSI INTERPOLASI DARI TABEL

from scipy import interpolate

def _func_interp_bivariate(df):
    "Membuat fungsi dari tabel untuk interpolasi bilinear"
    table = df[df.columns.sort_values()].sort_index().copy()

    x = table.index
    y = table.columns
    z = table.to_numpy()

    # penggunaan kx=1, ky=1 untuk interpolasi linear antara 2 titik
    # tidak menggunakan (cubic) spline interpolation
    return interpolate.RectBivariateSpline(x, y, z, kx=1, ky=1)

def _as_value(x):
    x = np.around(x, 4)
    return x.flatten() if x.size > 1 else x.item()
In [7]:
def find_K(prob, skew_log, source='scipy'):
    "Mencari nilai K berdasarkan probabilitas dan skewness logaritmik"
    prob = np.array(prob)
    
    if source.lower() == 'scipy':
        #ref: https://github.com/hidrokit/hidrokit/discussions/156
        if skew_log >= 0: 
            return np.around(stats.pearson3.ppf(1-prob, skew_log), 4)
        else:
            return np.around(stats.pearson3.ppf(prob, skew_log), 4)
    elif source.lower() == 'soetopo':
        func_pearson3_st = _func_interp_bivariate(t_pearson3_st)
        return _as_value(func_pearson3_st(skew_log, prob, grid=False))
    elif source.lower() == 'soewarno':
        func_pearson3_sw = _func_interp_bivariate(t_pearson3_sw)
        return _as_value(func_pearson3_sw(skew_log, prob, grid=False))
    elif source.lower() == 'limantara':
        func_pearson3_lm = _func_interp_bivariate(t_pearson3_lm)
        return _as_value(func_pearson3_lm(skew_log, prob, grid=False))

def calc_x_lp3(x, return_period=[5], source='scipy', show_stat=False):
    "Menghitung besar X dengan kala ulang tertentu"
    y = np.log10(x)
    y_mean = np.mean(y)
    y_std = np.std(y, ddof=1)
    y_skew = stats.skew(y, bias=False)

    prob = 1 / np.array(return_period)
    k = find_K(prob, y_skew, source=source)

    if show_stat:
        print(f'y_mean = {y_mean:.5f}')
        print(f'y_std = {y_std:.5f}')
        print(f'y_skew = {y_skew:.5f}')
        print(f'k = {k}')
    
    val_y = y_mean + k * y_std
    val_x = np.power(10, val_y)
    return val_x

def freq_logpearson3(
    df, col=None, 
    return_period=[2, 5, 10, 20, 25, 50, 100], source='scipy', show_stat=False,
    col_name='Log Pearson III', index_name='Kala Ulang'
    ):

    col = df.columns[0] if col is None else col

    x = df[col].copy()

    arr = calc_x_lp3(
        x, return_period=return_period, show_stat=show_stat,
        source=source
    )

    result = pd.DataFrame(
        data=arr, index=return_period, columns=[col_name]
    )

    result.index.name = index_name
    return result
In [8]:
dict_table_source = {
    'soewarno': t_pearson3_sw,
    'soetopo': t_pearson3_st,
    'limantara': t_pearson3_lm
}

def _find_prob_in_table(k, skew_log, table):
    func_table = _func_interp_bivariate(table)
    y = table.columns
    x = func_table(skew_log, y, grid=False)
    func_prob = interpolate.interp1d(x, y, kind='linear')
    return _as_value(func_prob(k))

def _calc_prob_in_table(k, skew_log, source='soewarno'):
    if source.lower() in dict_table_source.keys():
        return 1 - _find_prob_in_table(
            k, skew_log, dict_table_source[source.lower()]
        )
    
def calc_prob(k, skew_log, source='scipy'):
    if source.lower() == 'scipy':
        # TODO: CHECK/TEST THIS CODE
        if skew_log >= 0:
            return stats.pearson3.cdf(k, skew_log)
        else:
            return stats.pearson3.cdf(k, skew_log)
    if source.lower() in dict_table_source.keys():
        return _calc_prob_in_table(k, skew_log, source.lower())

FUNGSI

Fungsi find_K(prob, skew_log, ...)

Function: find_K(prob, skew_log, source='scipy')

Fungsi find_K(...) digunakan untuk mencari nilai $K$ (frequency factor) dari berbagai sumber berdasarkan besar probabilitas dan kemencengan (skew) logaritmik data.

  • Argumen Posisi:
    • prob: probabilitas $\left( \left( 0, 1 \right) \in \mathbb{R} \right)$
    • skew_log: skewness logaritmik
  • Argumen Opsional:
    • source: sumber nilai K, 'scipy' (default). Sumber yang dapat digunakan antara lain: Soewarno ('soewarno'), Soetopo ('soetopo'), Limantara ('limantara').
In [9]:
find_K(0.02, -2) # menggunakan nilai dari scipy
Out[9]:
0.9798
In [10]:
find_K(0.02, -2, source='soewarno') # menggunakan tabel dari buku Soewarno
Out[10]:
0.98
In [11]:
# perbandingan antara nilai tabel dan fungsi scipy

source_test = ['scipy', 'soewarno', 'soetopo', 'limantara']

prob = 1 / np.array([2, 5, 10, 20, 25, 50, 100]) # [0.5 , 0.2 , 0.1 , 0.05, 0.04, 0.02, 0.01]
print(f'prob = {prob}')
skew_log = -2.5
print(f'Nilai dari tabel (skew={skew_log})')
for _source in source_test:
    print(f'K_{_source:<12}=', find_K(prob, skew_log, source=_source))

skew_log = 1.75
print(f'Nilai hasil interpolasi (skew={skew_log})')
for _source in source_test:
    print(f'K_{_source:<12}=', find_K(prob, skew_log, source=_source))
prob = [0.5  0.2  0.1  0.05 0.04 0.02 0.01]
Nilai dari tabel (skew=-2.5)
K_scipy       = [0.3599 0.7107 0.7706 0.7901 0.7931 0.7977 0.7992]
K_soewarno    = [0.36   0.711  0.771  0.7893 0.793  0.798  0.799 ]
K_soetopo     = [0.36   0.711  0.771  0.7893 0.793  0.798  0.799 ]
K_limantara   = [0.36   0.711  0.771  0.7893 0.793  0.796  0.799 ]
Nilai hasil interpolasi (skew=1.75)
K_scipy       = [-0.2748  0.6515  1.3208  1.9769  2.1862  2.8318  3.472 ]
K_soewarno    = [-0.275   0.651   1.3207  2.0414  2.1855  2.831   3.4712]
K_soetopo     = [-0.275   0.6515  1.321   2.0418  2.186   2.8315  3.4715]
K_limantara   = [-0.275   0.6515  1.321   2.0418  2.186   2.8315  3.4715]

Fungsi calc_x_lp3(x, ...)

Function: calc_x_lp3(x, return_period=[5], source='scipy', show_stat=False)

Fungsi calc_x_lp3(...) digunakan untuk mencari besar $X$ berdasarkan kala ulang (return period), yang hasilnya dalam bentuk numpy.array.

  • Argumen Posisi:
    • x: array.
  • Argumen Opsional:
    • return_period: Kala Ulang (Tahun), [5] (default).
    • source: sumber nilai K, 'scipy' (default). Sumber yang dapat digunakan antara lain: Soewarno ('soewarno'), Soetopo ('soetopo'), Limantara ('limantara').
    • show_stat: Menampilkan parameter statistik, False (default)
In [12]:
calc_x_lp3(data.hujan)
Out[12]:
array([114.75316488])
In [13]:
calc_x_lp3(data.hujan, show_stat=True)
y_mean = 1.96408
y_std = 0.11381
y_skew = -1.25646
k = [0.8408]
Out[13]:
array([114.75316488])
In [14]:
calc_x_lp3(data.hujan, return_period=[5, 10, 15, 20, 21], show_stat=True)
y_mean = 1.96408
y_std = 0.11381
y_skew = -1.25646
k = [0.8408 1.0738 1.1678 1.2221 1.2304]
Out[14]:
array([114.75316488, 121.97804293, 125.02000217, 126.81166723,
       127.08778456])

Fungsi freq_logpearson3(df, ...)

Function: freq_logpearson3(df, col=None, return_period=[2, 5, 10, 20, 25, 50, 100], source='scipy', show_stat=False, col_name='Log Pearson III')

Fungsi freq_logpearson3(...) merupakan fungsi kembangan lebih lanjut dari calc_x_lp3(...) yang menerima input pandas.DataFrame dan memiliki luaran berupa pandas.DataFrame.

  • Argumen Posisi:
    • df: pandas.DataFrame.
  • Argumen Opsional:
    • col: nama kolom, None (default). Jika tidak diisi menggunakan kolom pertama dalam df sebagai data masukan.
    • return_period: Kala Ulang (Tahun), [2, 5, 10, 20, 25, 50, 100] (default).
    • source: sumber nilai K, 'scipy' (default). Sumber yang dapat digunakan antara lain: Soewarno ('soewarno'), Soetopo ('soetopo'), Limantara ('limantara').
    • show_stat: Menampilkan parameter statistik, False (default).
    • col_name: Nama kolom luaran, Log Pearson III (default).
In [15]:
freq_logpearson3(data)
Out[15]:
Log Pearson III
Kala Ulang
2 97.111291
5 114.753165
10 121.978043
20 126.811667
25 128.030417
50 131.065229
100 133.264423
In [16]:
freq_logpearson3(data, source='soewarno', col_name=f'LP3 (soewarno)')
Out[16]:
LP3 (soewarno)
Kala Ulang
2 97.103657
5 114.747151
10 121.962062
20 126.991242
25 128.020353
50 131.072098
100 133.281885
In [17]:
freq_logpearson3(data, 'hujan', source='limantara', col_name=f'LP3 (limantara)', show_stat=True)
y_mean = 1.96408
y_std = 0.11381
y_skew = -1.25646
k = [0.2035 0.8406 1.0736 1.2275 1.2583 1.3479 1.4117]
Out[17]:
LP3 (limantara)
Kala Ulang
2 97.103657
5 114.747151
10 121.971650
20 126.991242
25 128.020353
50 131.061795
100 133.271408
In [18]:
_res = []

for _s in ['scipy', 'soewarno', 'soetopo', 'limantara']:
    _res += [freq_logpearson3(data, 'hujan', source=_s, col_name=f'LP3 ({_s})')]

pd.concat(_res, axis=1)
Out[18]:
LP3 (scipy) LP3 (soewarno) LP3 (soetopo) LP3 (limantara)
Kala Ulang
2 97.111291 97.103657 97.103657 97.103657
5 114.753165 114.747151 114.747151 114.747151
10 121.978043 121.962062 121.971650 121.971650
20 126.811667 126.991242 126.991242 126.991242
25 128.030417 128.020353 128.020353 128.020353
50 131.065229 131.072098 131.061795 131.061795
100 133.264423 133.281885 133.271408 133.271408

Fungsi calc_prob(k, skew_log, ...)

Function: calc_prob(k, skew_log, source='scipy')

Fungsi calc_prob(...) digunakan untuk menghitung nilai probabilitas/peluang berdasarkan nilai $K$ (frequency factor).

  • Argumen Posisi:
    • k: nilai $K$ (frequency factor). Nilai $K$ diperoleh menggunakan persamaan $K = \frac{y-\bar{y}}{s_y}$ dengan $y = log(x)$.
    • skew_log: nilai skewness logaritmik.
  • Argumen Opsional:
    • source: sumber nilai probabilitas. 'scipy' (default). Sumber yang dapat digunakan antara lain: Soewarno ('soewarno'), Soetopo ('soetopo'), Limantara ('limantara'). Ketiga sumber lain menggunakan tabel, sehingga memiliki keterbatasan dalam memberi nilai probabilitas.

Probabilitas/Peluang yang dikeluarkan oleh fungsi calc_prob(...) adalah $P'(X<)$ sehingga jika dilakukan perhitungan $K$ kembali, nilai masukan harus menggunakan formula $P(X) = 1-P'(X<)$.

In [19]:
_data = np.array([59.000, 67.750, 50.000, 72.500, 29.000, 30.550, 36.100, 64.500, 34.500, 125.500])
_y = np.sort(np.log10(_data))[::-1] #urut besar ke kecil
_y_mean = np.mean(_y)
_y_std = np.std(_y, ddof=1)
_y_skew = stats.skew(_y, bias=False)

_k = (_y - _y_mean) / _y_std
_k
Out[19]:
array([ 1.91913743,  0.73872755,  0.59295708,  0.48720527,  0.29547275,
       -0.06058346, -0.76129719, -0.85881912, -1.12039465, -1.23240566])
In [20]:
calc_prob(_k, _y_skew, source='limantara')
Out[20]:
array([0.9607, 0.7761, 0.727 , 0.6913, 0.6266, 0.5065, 0.2365, 0.1992,
       0.1271, 0.0975])
In [21]:
calc_prob(_k, _y_skew, source='soetopo')
Out[21]:
array([0.9607, 0.7761, 0.727 , 0.6913, 0.6267, 0.5066, 0.2368, 0.1994,
       0.1272, 0.0975])
In [22]:
calc_prob(_k, _y_skew, source='scipy')
Out[22]:
array([0.96095003, 0.78178546, 0.7409749 , 0.70849576, 0.64381381,
       0.50782301, 0.23212832, 0.19918497, 0.12334252, 0.09705203])

Changelog

- 20220613 - 1.2.0 / v0.4.0 - ubah untuk perhitungan calc_prob (CDF) sehingga lebih sederhana dan tidak menggunakan negasi (tidak dipengaruhi skewness). Perlu diteliti lebih lanjut. 
- 20220323 - 1.1.0 - tambah argumen index_name="Kala Ulang" pada fungsi freq_gumbel() untuk penamaan index
- 20220315 - 1.0.3 - Tambah fungsi calc_prob(...)
- 20220310 - 1.0.2 - Fix show_stat default typo
- 20220309 - 1.0.1 - Typo
- 20220309 - 1.0.0 - Initial

Copyright © 2022 Taruma Sakti Megariansyah

Source code in this notebook is licensed under a MIT License. Data in this notebook is licensed under a Creative Common Attribution 4.0 International.