library(samurais)
## Warning: package 'samurais' was built under R version 3.6.1

Package‘samurais’ in R (2019 best packages for time series forecasting)

Provides a variety of original and flexible user-friendly statistical latent variable models and unsupervised learning algorithms to segment,estimate ,forecast and represent time-series HPI data (univariate or multivariate), and more generally, longitudinal data, which include regime changes.

‘samurais package in R’ is built upon the following packages, each of them is an autonomous time-series segmentation approach: Regression with Hidden Logistic Process (‘RHLP’), Hidden Markov Model Regression (‘HMMR’), Multivariate ‘RHLP’ (‘MRHLP’), Multivariate ‘HMMR’ (‘MHMMR’), Piece-Wise regression (‘PWR’).

library(readr)
## Warning: package 'readr' was built under R version 3.6.1
alaska<- read_csv("~/RPy Analytics/Dura capital/hpi_data.csv")
## Parsed with column specification:
## cols(
##   .default = col_double()
## )
## See spec(...) for full column specifications.
#View(hpi_data)
## univariate scenario (Alaska)
data(univtoydataset) 
#univtoydataset
hmmr <- emHMMR(alaska$year, alaska$AL, K = 10, p = 1, verbose = TRUE)
## EM: Iteration : 1 || log-likelihood : -766.832644939466
## EM: Iteration : 2 || log-likelihood : -535.504623683392
## EM: Iteration : 3 || log-likelihood : -502.119297097182
## EM: Iteration : 4 || log-likelihood : -483.580586739892
## EM: Iteration : 5 || log-likelihood : -480.50231409789
## EM: Iteration : 6 || log-likelihood : -478.917779395797
## EM: Iteration : 7 || log-likelihood : -476.329126841241
## EM: Iteration : 8 || log-likelihood : -473.104679006341
## EM: Iteration : 9 || log-likelihood : -470.156955966773
## EM: Iteration : 10 || log-likelihood : -468.849840400063
## EM: Iteration : 11 || log-likelihood : -468.712247491231
## EM: Iteration : 12 || log-likelihood : -468.688383897746
## EM: Iteration : 13 || log-likelihood : -468.677059829944
## EM: Iteration : 14 || log-likelihood : -468.671256037552
## EM: Iteration : 15 || log-likelihood : -468.668300155241
## EM: Iteration : 16 || log-likelihood : -468.666805535232
## EM: Iteration : 17 || log-likelihood : -468.666046888169
## EM: Iteration : 18 || log-likelihood : -468.665656224047
hmmr$summary()
## ---------------------
## Fitted HMMR model
## ---------------------
## 
## HMMR model with K = 10 components:
## 
##  log-likelihood  nu       AIC       BIC
##       -468.6657 129 -597.6657 -802.5273
## 
## Clustering table (Number of observations in each regimes):
## 
##  1  2  3  4  5  6  7  8  9 10 
## 18 12 17 18 26 11 20 16 18 21 
## 
## Regression coefficients:
## 
##       Beta(K = 1)  Beta(K = 2)   Beta(K = 3)  Beta(K = 4)  Beta(K = 5)
## 1   -13766.808877 -3842.150469 -10373.569410 -5274.597315 -12698.28782
## X^1      7.006813     1.994201      5.289747     2.724585      6.45406
##     Beta(K = 6)   Beta(K = 7)   Beta(K = 8) Beta(K = 9) Beta(K = 10)
## 1   -9695.07113 -17512.503369 -16465.752633 9806.980524 -20591.06466
## X^1     4.95302      8.862867      8.350838   -4.736351     10.36144
## 
## Variances:
## 
##  Sigma2(K = 1) Sigma2(K = 2) Sigma2(K = 3) Sigma2(K = 4) Sigma2(K = 5)
##        6.91034      12.30262      5.641697      1.540016      4.122403
##  Sigma2(K = 6) Sigma2(K = 7) Sigma2(K = 8) Sigma2(K = 9) Sigma2(K = 10)
##       3.165992      8.684245       32.6222      7.544492       15.72537
hmmr$plot()

### multvariate scenario
mhmmr <- emMHMMR(alaska$year,alaska[,c("AK",    "AL",   "AR",   "AZ",   "CA")], K = 10, p = 1, verbose = TRUE)
## EM: Iteration : 1 || log-likelihood : -4392.39312424105
## EM: Iteration : 2 || log-likelihood : -2280.49809640766
## EM: Iteration : 3 || log-likelihood : -2275.4110009901
## EM: Iteration : 4 || log-likelihood : -2272.49913062106
## EM: Iteration : 5 || log-likelihood : -2270.13551918563
## EM: Iteration : 6 || log-likelihood : -2269.81801187087
## EM: Iteration : 7 || log-likelihood : -2269.21646549831
## EM: Iteration : 8 || log-likelihood : -2268.76855045316
## EM: Iteration : 9 || log-likelihood : -2268.66239434338
## EM: Iteration : 10 || log-likelihood : -2268.5173684659
## EM: Iteration : 11 || log-likelihood : -2267.91195069571
## EM: Iteration : 12 || log-likelihood : -2266.01773290908
## EM: Iteration : 13 || log-likelihood : -2265.54002211479
## EM: Iteration : 14 || log-likelihood : -2265.53997367785
mhmmr$summary()
## ----------------------
## Fitted MHMMR model
## ----------------------
## 
## MHMMR model with K = 10 regimes
## 
##  log-likelihood  nu      AIC       BIC
##        -2265.54 349 -2614.54 -3168.778
## 
## Clustering table:
##  1  2  3  4  5  6  7  8  9 10 
## 17 17 17 15 19 19 15 16 18 24 
## 
## 
## ------------------
## Regime 1 (K = 1):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   1.933309e-05 2.124384e-05 1.916561e-05 1.700387e-05 0.0000155451
## X^1 3.821441e-02 4.199128e-02 3.788337e-02 3.361039e-02 0.0307269546
## 
## Covariance matrix:
##                                                  
##   61.88105  64.87395  65.82106  60.31517 104.9246
##   64.87395  76.94748  76.67896  73.18855 121.0847
##   65.82106  76.67896  79.28542  72.55781 123.0497
##   60.31517  73.18855  72.55781  78.18367 114.7617
##  104.92456 121.08471 123.04968 114.76175 198.0223
## ------------------
## Regime 2 (K = 2):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   2.827447e-05 2.754261e-05 2.708866e-05 2.725869e-05 2.777031e-05
## X^1 5.600808e-02 5.455835e-02 5.365913e-02 5.399596e-02 5.500939e-02
## 
## Covariance matrix:
##                                                 
##  480.32550 52.34772 55.92119 100.51901 140.47317
##   52.34772 20.95463 20.26699  30.00160  32.34991
##   55.92119 20.26699 29.17683  37.77143  40.12183
##  100.51901 30.00160 37.77143  61.44997  67.31065
##  140.47317 32.34991 40.12183  67.31065  89.31451
## ------------------
## Regime 3 (K = 3):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   3.558194e-05 3.246503e-05 3.149696e-05 3.230926e-05 3.319651e-05
## X^1 7.063436e-02 6.444691e-02 6.252519e-02 6.413771e-02 6.589899e-02
## 
## Covariance matrix:
##                                                   
##   88.99223 -54.86067 -28.09875 -48.39661 -90.23697
##  -54.86067  55.81858  28.49121  52.40100  76.35080
##  -28.09875  28.49121  16.69622  28.30399  40.00428
##  -48.39661  52.40100  28.30399  55.28000  72.54715
##  -90.23697  76.35080  40.00428  72.54715 112.89830
## ------------------
## Regime 4 (K = 4):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   2.900766e-05 0.0000366562 3.235088e-05 3.385254e-05 5.177285e-05
## X^1 5.770013e-02 0.0729141096 6.435024e-02 6.733725e-02 1.029832e-01
## 
## Covariance matrix:
##                                                       
##  111.026581  1.6044571  3.8107743  2.9773311 -41.89270
##    1.604457 10.2219346  3.9745888  0.5408306  71.31540
##    3.810774  3.9745888  1.8145033  0.4108295  26.25097
##    2.977331  0.5408306  0.4108295  1.1041686  -6.91685
##  -41.892696 71.3154049 26.2509674 -6.9168503 678.07022
## ------------------
## Regime 5 (K = 5):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   3.501514e-05 4.201636e-05 3.672467e-05 3.731638e-05 5.374788e-05
## X^1 6.979823e-02 8.375429e-02 7.320596e-02 7.438546e-02 1.071396e-01
## 
## Covariance matrix:
##                                                   
##   46.60667  60.52686  63.97621  59.44579 -65.48396
##   60.52686  81.17453  84.81540  77.71562 -86.10664
##   63.97621  84.81540  89.75605  82.33438 -92.05009
##   59.44579  77.71562  82.33438  77.44808 -81.83248
##  -65.48396 -86.10664 -92.05009 -81.83248 105.93579
## ------------------
## Regime 6 (K = 6):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   4.115695e-05 5.006855e-05 4.352742e-05 4.632426e-05 5.786691e-05
## X^1 8.223606e-02 1.000424e-01 8.697252e-02 9.256091e-02 1.156244e-01
## 
## Covariance matrix:
##                                                  
##   38.08685  57.74382  41.41296  77.19636 155.3075
##   57.74382  90.62697  64.68069 118.33223 236.1848
##   41.41296  64.68069  46.46492  85.36845 171.5960
##   77.19636 118.33223  85.36845 162.03042 331.5313
##  155.30748 236.18483 171.59596 331.53126 692.8297
## ------------------
## Regime 7 (K = 7):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   4.883203e-05 5.832711e-05 5.059906e-05 5.909917e-05 9.155577e-05
## X^1 9.778129e-02 1.167942e-01 1.013196e-01 1.183402e-01 1.833313e-01
## 
## Covariance matrix:
##                                                  
##  182.3675 113.07680 121.03686  245.4053  748.5601
##  113.0768  75.15479  78.40509  153.5705  469.9886
##  121.0369  78.40509  82.97036  164.3159  502.3604
##  245.4053 153.57046 164.31589  333.5739 1017.5771
##  748.5601 469.98862 502.36041 1017.5771 3107.8847
## ------------------
## Regime 8 (K = 8):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   6.638267e-05 7.103192e-05 6.141522e-05 9.540986e-05 0.0001449247
## X^1 1.331803e-01 1.425078e-01 1.232143e-01 1.914161e-01 0.2907553427
## 
## Covariance matrix:
##                                                     
##  312.7745 295.938209 188.82385  564.0075  184.034515
##  295.9382 292.286889 180.60933  451.1846    1.650569
##  188.8239 180.609325 114.63278  333.2994   94.926965
##  564.0075 451.184554 333.29943 1729.5584 1799.790899
##  184.0345   1.650569  94.92697 1799.7909 3179.904983
## ------------------
## Regime 9 (K = 9):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4) Beta(d = 5)
## 1   6.933623e-05 7.073846e-05 6.096503e-05 6.532265e-05 9.97045e-05
## X^1 1.394006e-01 1.422198e-01 1.225703e-01 1.313313e-01 2.00456e-01
## 
## Covariance matrix:
##                                                      
##   6.859675  -5.77749   1.560809   3.991844   3.031998
##  -5.777490 104.78908  46.212724 287.948580 170.547741
##   1.560809  46.21272  23.598190 142.584729  84.044967
##   3.991844 287.94858 142.584729 907.838460 534.848652
##   3.031998 170.54774  84.044967 534.848652 336.633852
## ------------------
## Regime 10 (K = 10):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   7.623422e-05 7.285109e-05 6.514573e-05 0.0000847447 0.0001370561
## X^1 1.536693e-01 1.468497e-01 1.313176e-01 0.1708242987 0.2762711039
## 
## Covariance matrix:
##                                                  
##  139.9193  192.2679  162.2606  502.0126  830.9939
##  192.2679  289.7827  243.1845  739.3775 1193.0375
##  162.2606  243.1845  206.4592  625.3424 1010.7015
##  502.0126  739.3775  625.3424 1920.5276 3124.5251
##  830.9939 1193.0375 1010.7015 3124.5251 5127.9929
mhmmr$plot()

“AK”, “AL”, “AR”, “AZ”, “CA”, “CO”,“CT”,“DE”,“FL”,“GA”,“HI”,“IA”,“ID”,“IL”,“IN”,“KS”,“KY”,“LA”,“MA”,“MD”,“ME”,“MI”,“MN”, “MO”,“MS”,“MT”,“NC”,“ND”,“NE”,“NH”,“NJ”,“NM”,“NV”,“NY”,“OH”,“OK”,“OR”,“PA”,“RI”,“SC”,“SD”, “TN”,“TX”,“UT”,“VA”,“VT”,“WA”,“WI”,“WV”,“WY”,“DC”,“US”

##emMRHLP implements the maximum-likelihood parameter estimation of the MRHLP model by the Expectation-Maximization (EM) algorithm.

mrhlp <- emMRHLP(alaska$year,alaska[,c("AK",    "AL",   "AR",   "AZ",   "CA")], K = 10, p = 1, verbose = TRUE)
## EM: Iteration : 1 || log-likelihood : -4741.31572216474
## EM: Iteration : 2 || log-likelihood : -2003.29238249833
## EM: Iteration : 3 || log-likelihood : -2003.13543888782
## EM: Iteration : 4 || log-likelihood : -2002.93931844023
## EM: Iteration : 5 || log-likelihood : -2002.58710746702
## EM: Iteration : 6 || log-likelihood : -1997.17812068174
## EM: Iteration : 7 || log-likelihood : -1994.88312737565
## EM: Iteration : 8 || log-likelihood : -1994.03349478655
## EM: Iteration : 9 || log-likelihood : -1993.96673787051
## EM: Iteration : 10 || log-likelihood : -1993.96494110611
mrhlp$summary()
## ----------------------
## Fitted MRHLP model
## ----------------------
## 
## MRHLP model with K = 10 regimes
## 
##  log-likelihood  nu       AIC       BIC       ICL
##       -1993.965 268 -2261.965 -2687.569 -2816.147
## 
## Clustering table:
##  1  2  3  4  5  6  7  8  9 10 
## 16 20 16 16 16 20 16 16 16 25 
## 
## 
## ------------------
## Regime 1 (K = 1):
## 
## Regression coefficients:
## 
##       Beta(d = 1)   Beta(d = 2)   Beta(d = 3)   Beta(d = 4)  Beta(d = 5)
## 1   -11900.904154 -13504.404091 -13859.981518 -12877.163165 -21967.01632
## X^1      6.059007      6.874019      7.049803      6.548325     11.14407
## 
## Covariance matrix:
##                                                
##  7.306377 3.170526 2.534520  1.879064  4.303645
##  3.170526 7.164631 5.121944  7.169621  7.258622
##  2.534520 5.121944 5.892427  4.877967  6.304251
##  1.879064 7.169621 4.877967 15.834046  6.973094
##  4.303645 7.258622 6.304251  6.973094 12.375397
## ------------------
## Regime 2 (K = 2):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2) Beta(d = 3)   Beta(d = 4)   Beta(d = 5)
## 1   -27616.16229 -4873.247560 -6404.82360 -10124.997781 -12818.871824
## X^1     13.99742     2.514723     3.28702      5.165413      6.526377
## 
## Covariance matrix:
##                                                      
##  162.521049 -4.489663 -19.240596 -17.995227 -9.471834
##   -4.489663 10.791472   6.725710   8.674966  5.327508
##  -19.240596  6.725710  11.160205   9.393123  4.167246
##  -17.995227  8.674966   9.393123  16.832009 10.790537
##   -9.471834  5.327508   4.167246  10.790537 17.826252
## ------------------
## Regime 3 (K = 3):
## 
## Regression coefficients:
## 
##      Beta(d = 1)   Beta(d = 2)  Beta(d = 3)   Beta(d = 4)   Beta(d = 5)
## 1   11509.690199 -11059.552668 -5896.513982 -11190.784916 -15683.939259
## X^1    -5.727352      5.635678     3.032884      5.701477      7.966657
## 
## Covariance matrix:
##                                                          
##   34.11026952 -2.113539 0.02422109  4.9800251 -15.4398628
##   -2.11353864  5.096288 1.45041317  1.0858916   4.4376027
##    0.02422109  1.450413 2.28054863  0.9473013   1.6650559
##    4.98002509  1.085892 0.94730133  3.3619984  -0.2113683
##  -15.43986282  4.437603 1.66505591 -0.2113683  10.9351270
## ------------------
## Regime 4 (K = 4):
## 
## Regression coefficients:
## 
##       Beta(d = 1)  Beta(d = 2)  Beta(d = 3) Beta(d = 4) Beta(d = 5)
## 1   -1593.9182685 -5055.699778 -2087.281190 98.16118069 -40493.6241
## X^1     0.8590119     2.614572     1.113692  0.01798848     20.4604
## 
## Covariance matrix:
##                                                        
##  110.169732 -1.0886254  2.6976560  3.0277559 -63.409161
##   -1.088625  1.7149262  0.4617160  0.7051391   3.205330
##    2.697656  0.4617160  0.3637622  0.4785143  -1.870394
##    3.027756  0.7051391  0.4785143  1.1009042  -5.597524
##  -63.409161  3.2053302 -1.8703939 -5.5975241 132.643579
## ------------------
## Regime 5 (K = 5):
## 
## Regression coefficients:
## 
##      Beta(d = 1)   Beta(d = 2)   Beta(d = 3)   Beta(d = 4)  Beta(d = 5)
## 1   -9209.856535 -12265.778674 -12865.753944 -11750.906964 13525.258262
## X^1     4.690035      6.237032      6.527469      5.969372    -6.677972
## 
## Covariance matrix:
##                                                     
##   3.286104  2.846715  3.475051  4.1814685 -1.8874939
##   2.846715  4.372768  4.256318  4.1316165 -1.4221758
##   3.475051  4.256318  5.256036  5.1509030 -3.2229164
##   4.181469  4.131617  5.150903  6.9468166 -0.6972664
##  -1.887494 -1.422176 -3.222916 -0.6972664 12.5610274
## ------------------
## Regime 6 (K = 6):
## 
## Regression coefficients:
## 
##      Beta(d = 1)   Beta(d = 2)  Beta(d = 3)   Beta(d = 4)  Beta(d = 5)
## 1   -8554.177978 -13171.414066 -9367.802442 -17370.205638 -33987.53966
## X^1     4.363342      6.691917     4.775236      8.785709     17.12488
## 
## Covariance matrix:
##                                                
##  2.650534 2.810083 1.891754  3.086018  5.535062
##  2.810083 5.616713 3.476814  3.392068  3.508881
##  1.891754 3.476814 2.378730  2.526466  3.562539
##  3.086018 3.392068 2.526466  6.524568 15.711286
##  5.535062 3.508881 3.562539 15.711286 48.063374
## ------------------
## Regime 7 (K = 7):
## 
## Regression coefficients:
## 
##     Beta(d = 1)   Beta(d = 2)   Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   -21998.3367 -15434.212546 -15694.650635 -29652.91384 -89519.92513
## X^1     11.0839      7.824645      7.939285     14.92723     44.89043
## 
## Covariance matrix:
##                                                    
##   25.635793  8.181669 11.850899  33.25418 100.87891
##    8.181669  5.365530  5.481485  11.33113  34.43508
##   11.850899  5.481485  6.920270  16.34983  49.95468
##   33.254178 11.331125 16.349834  46.23862 140.62418
##  100.878911 34.435080 49.954684 140.62418 433.13192
## ------------------
## Regime 8 (K = 8):
## 
## Regression coefficients:
## 
##     Beta(d = 1)  Beta(d = 2)   Beta(d = 3)  Beta(d = 4) Beta(d = 5)
## 1   -24713.6821 -25927.61905 -15043.966629 -17681.17066 39307.48881
## X^1     12.4513     13.06577      7.621596      9.00237   -19.30485
## 
## Covariance matrix:
##                                                   
##   74.05748  49.39038  42.81051  322.5282  432.8223
##   49.39038  37.17086  29.52319  198.0967  253.8118
##   42.81051  29.52319  25.34477  190.0608  255.3263
##  322.52819 198.09674 190.06084 1681.6615 2388.7545
##  432.82229 253.81183 255.32632 2388.7545 3492.1826
## ------------------
## Regime 9 (K = 9):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2) Beta(d = 3) Beta(d = 4)  Beta(d = 5)
## 1   -2429.760438 14061.497486 5040.268423 30358.48406 18372.196015
## X^1     1.347804    -6.851815   -2.384474   -14.96942    -8.938205
## 
## Covariance matrix:
##                                                  
##   4.898935  5.16892  5.479932  26.17117  16.27994
##   5.168920 16.74982 12.546168  72.24687  37.48382
##   5.479932 12.54617 10.653530  58.00309  31.57697
##  26.171169 72.24687 58.003092 352.01935 189.08025
##  16.279941 37.48382 31.576965 189.08025 122.65987
## ------------------
## Regime 10 (K = 10):
## 
## Regression coefficients:
## 
##       Beta(d = 1)   Beta(d = 2)   Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   -13077.228387 -18834.421123 -15859.339573 -49071.96938 -80435.94429
## X^1      6.641188      9.490471      7.999022     24.51507     40.17996
## 
## Covariance matrix:
##                                                   
##   9.292744  4.111581  3.826747  11.82554  27.54543
##   4.111581 18.767451 14.980264  33.30717  35.72764
##   3.826747 14.980264 14.302977  30.80712  36.20813
##  11.825536 33.307167 30.807118  81.05089 109.49598
##  27.545435 35.727636 36.208134 109.49598 186.17800
mrhlp$plot()

##univariate scenario

##emRHLP implements the maximum-likelihood parameter estimation of the RHLP model by the Expectation-Maximization (EM) algorithm
rhlp <- emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose = TRUE)
## EM: Iteration : 1 || log-likelihood : -777.043930149238
## EM: Iteration : 2 || log-likelihood : -568.750109903377
## EM: Iteration : 3 || log-likelihood : -541.495841533048
## EM: Iteration : 4 || log-likelihood : -524.398890765
## EM: Iteration : 5 || log-likelihood : -511.790453532024
## EM: Iteration : 6 || log-likelihood : -500.627209424373
## EM: Iteration : 7 || log-likelihood : -491.915032136676
## EM: Iteration : 8 || log-likelihood : -486.054970871925
## EM: Iteration : 9 || log-likelihood : -481.204341491301
## EM: Iteration : 10 || log-likelihood : -476.112648624409
## EM: Iteration : 11 || log-likelihood : -470.974242758004
## EM: Iteration : 12 || log-likelihood : -468.193275646603
## EM: Iteration : 13 || log-likelihood : -466.070291660498
## EM: Iteration : 14 || log-likelihood : -464.655538728975
## EM: Iteration : 15 || log-likelihood : -463.165834337997
## EM: Iteration : 16 || log-likelihood : -461.924585820573
## EM: Iteration : 17 || log-likelihood : -460.537682884993
## EM: Iteration : 18 || log-likelihood : -459.215737168283
## EM: Iteration : 19 || log-likelihood : -457.654385522431
## EM: Iteration : 20 || log-likelihood : -456.366865397772
## EM: Iteration : 21 || log-likelihood : -455.589333475138
## EM: Iteration : 22 || log-likelihood : -455.360012748036
## EM: Iteration : 23 || log-likelihood : -454.858292164273
## EM: Iteration : 24 || log-likelihood : -454.397213243937
## EM: Iteration : 25 || log-likelihood : -453.87336431414
## EM: Iteration : 26 || log-likelihood : -453.316540092209
## EM: Iteration : 27 || log-likelihood : -452.699082461163
## EM: Iteration : 28 || log-likelihood : -452.012586863723
## EM: Iteration : 29 || log-likelihood : -451.402481073369
## EM: Iteration : 30 || log-likelihood : -450.833099311636
## EM: Iteration : 31 || log-likelihood : -450.470665066122
## EM: Iteration : 32 || log-likelihood : -450.034435996985
## EM: Iteration : 33 || log-likelihood : -448.913580301846
## EM: Iteration : 34 || log-likelihood : -446.155958919772
## EM: Iteration : 35 || log-likelihood : -443.691397698499
## EM: Iteration : 36 || log-likelihood : -442.45784058546
## EM: Iteration : 37 || log-likelihood : -441.511300816952
## EM: Iteration : 38 || log-likelihood : -441.245680912435
## EM: Iteration : 39 || log-likelihood : -441.251634257156
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.245680912435to
## -441.251634257156 !
## EM: Iteration : 40 || log-likelihood : -441.250257193768
## EM: Iteration : 41 || log-likelihood : -441.261684696638
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.250257193768to
## -441.261684696638 !
## EM: Iteration : 42 || log-likelihood : -441.277704739432
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.261684696638to
## -441.277704739432 !
## EM: Iteration : 43 || log-likelihood : -441.294668572067
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.277704739432to
## -441.294668572067 !
## EM: Iteration : 44 || log-likelihood : -441.310263587088
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.294668572067to
## -441.310263587088 !
## EM: Iteration : 45 || log-likelihood : -441.326040402
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose = TRUE):
## EM log-likelihood is decreasing from -441.310263587088to -441.326040402 !
## EM: Iteration : 46 || log-likelihood : -441.335404635877
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose = TRUE):
## EM log-likelihood is decreasing from -441.326040402to -441.335404635877 !
## EM: Iteration : 47 || log-likelihood : -441.342271119417
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.335404635877to
## -441.342271119417 !
## EM: Iteration : 48 || log-likelihood : -441.346540861214
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.342271119417to
## -441.346540861214 !
## EM: Iteration : 49 || log-likelihood : -441.350463827192
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.346540861214to
## -441.350463827192 !
## EM: Iteration : 50 || log-likelihood : -441.352498283004
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.350463827192to
## -441.352498283004 !
## EM: Iteration : 51 || log-likelihood : -441.355032491954
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.352498283004to
## -441.355032491954 !
## EM: Iteration : 52 || log-likelihood : -441.356874356842
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.355032491954to
## -441.356874356842 !
## EM: Iteration : 53 || log-likelihood : -441.359107163707
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.356874356842to
## -441.359107163707 !
## EM: Iteration : 54 || log-likelihood : -441.364998758417
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.359107163707to
## -441.364998758417 !
## EM: Iteration : 55 || log-likelihood : -441.374307523967
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.364998758417to
## -441.374307523967 !
## EM: Iteration : 56 || log-likelihood : -441.359615504697
## EM: Iteration : 57 || log-likelihood : -441.374552023265
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.359615504697to
## -441.374552023265 !
## EM: Iteration : 58 || log-likelihood : -441.392086987475
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.374552023265to
## -441.392086987475 !
## EM: Iteration : 59 || log-likelihood : -441.409935140254
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.392086987475to
## -441.409935140254 !
## EM: Iteration : 60 || log-likelihood : -441.425826244974
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.409935140254to
## -441.425826244974 !
## EM: Iteration : 61 || log-likelihood : -441.442106738065
## Warning in emRHLP(alaska$year, alaska$AL, K = 10, p = 2, verbose
## = TRUE): EM log-likelihood is decreasing from -441.425826244974to
## -441.442106738065 !
rhlp$summary()
## ---------------------
## Fitted RHLP model
## ---------------------
## 
## RHLP model with K = 10 components:
## 
##  log-likelihood nu       AIC       BIC       ICL
##       -441.4421 58 -499.4421 -591.5504 -710.5483
## 
## Clustering table (Number of observations in each regimes):
## 
##  1  2  3  4  5  7  8  9 10 
##  8  8 48 48  4 12  8 20 21 
## 
## Regression coefficients:
## 
##      Beta(K = 1)   Beta(K = 2)   Beta(K = 3)   Beta(K = 4)   Beta(K = 5)
## 1   82.032127362 -94.767726200 -1.981986e+04  1.068310e+05  2.730778e+03
## X^1 -4.974627831 -16.192289931  1.563457e+01 -1.136945e+02 -8.080012e+00
## X^2  0.002516369   0.008166201 -2.813833e-03  3.019285e-02  3.409541e-03
##       Beta(K = 6)   Beta(K = 7)   Beta(K = 8)   Beta(K = 9)  Beta(K = 10)
## 1   -3.607231e+02 -2.959367e+03 -2.564165e+03  1.454670e+04  1.704940e+04
## X^1 -6.755822e+00 -5.220290e+00  6.955369e-02 -9.434273e+00 -2.702613e+01
## X^2  3.521719e-03  3.403794e-03  6.763459e-04  1.164110e-03  9.284066e-03
## 
## Variances:
## 
##  Sigma2(K = 1) Sigma2(K = 2) Sigma2(K = 3) Sigma2(K = 4) Sigma2(K = 5)
##       1.469476      19.00826      9.724914      4.972181      1.847543
##  Sigma2(K = 6) Sigma2(K = 7) Sigma2(K = 8) Sigma2(K = 9) Sigma2(K = 10)
##     0.01776632      2.138432      4.724555      7.440869       15.48246
rhlp$plot()

## univariate Alasak

###fitPWRFisherisusedtofitaPiecewiseRegression(PWR)modelbymaximum-likelihoodviaanoptimized dynamic programming algorithm. The estimation performed by the dynamic programming algorithm provides an optimal segmentation of the time series.

pwr <- fitPWRFisher(alaska$year, alaska$AL, K = 5, p = 2)
pwr$summary()
## --------------------
## Fitted PWR model
## --------------------
## 
## PWR model with K = 5 components:
## 
## Clustering table (Number of observations in each regimes):
## 
##  1  2  3  4  5 
## 60 58 14 24 21 
## 
## Regression coefficients:
## 
##       Beta(K = 1)   Beta(K = 2)   Beta(K = 3)   Beta(K = 4)   Beta(K = 5)
## 1   -3.409539e+05  3.631015e+05 -1.111570e+06  3.747469e+06  2.852111e+06
## X^1  3.390493e+02 -3.706222e+02  1.091943e+03 -3.721325e+03 -2.839329e+03
## X^2 -8.424163e-02  9.458924e-02 -2.680354e-01  9.239084e-01  7.067152e-01
## 
## Variances:
## 
##  Sigma2(K = 1) Sigma2(K = 2) Sigma2(K = 3) Sigma2(K = 4) Sigma2(K = 5)
##       11.78053      5.271198      25.74004       10.0187      13.19492
pwr$plot()

## Alaska

###hmmProcess calculates the probability distribution of a random process following a Markov chain of HPIs overtime

hmmr <- emHMMR(alaska$year, alaska$AL, K = 5, p = 1)
# hmmr is a ModelHMMR object. It contains some methods such as 'summary' and 'plot'
hmmr$summary() 
## ---------------------
## Fitted HMMR model
## ---------------------
## 
## HMMR model with K = 5 components:
## 
##  log-likelihood nu       AIC       BIC
##       -549.5994 39 -588.5994 -650.5343
## 
## Clustering table (Number of observations in each regimes):
## 
##  1  2  3  4  5 
## 47 19 37 19 55 
## 
## Regression coefficients:
## 
##       Beta(K = 1)  Beta(K = 2)   Beta(K = 3)   Beta(K = 4)  Beta(K = 5)
## 1   -10207.822536 -5298.951730 -13540.762973 -17381.873931 -2410.424680
## X^1      5.206819     2.736838      6.876756      8.797582     1.343093
## 
## Variances:
## 
##  Sigma2(K = 1) Sigma2(K = 2) Sigma2(K = 3) Sigma2(K = 4) Sigma2(K = 5)
##       15.70333      1.554445      5.447689      8.792998      190.7014
hmmr$plot()

# hmmr has also two fields, stat and param which are reference classes as well
# Log-likelihood: 
hmmr$stat$loglik
## [1] -549.5994
# Parameters of the polynomial regressions: 
hmmr$param$beta
##               [,1]         [,2]          [,3]          [,4]         [,5]
## [1,] -10207.822536 -5298.951730 -13540.762973 -17381.873931 -2410.424680
## [2,]      5.206819     2.736838      6.876756      8.797582     1.343093
### multivariate 

####ModelMHMMR represents an estimated MHMMR model
mhmmr <- emMHMMR(alaska$year,alaska[,c("AK",    "AL",   "AR",   "AZ",   "CA")], K = 10, p = 1, verbose = TRUE)
## EM: Iteration : 1 || log-likelihood : -4392.39312424105
## EM: Iteration : 2 || log-likelihood : -2280.49809640766
## EM: Iteration : 3 || log-likelihood : -2275.4110009901
## EM: Iteration : 4 || log-likelihood : -2272.49913062106
## EM: Iteration : 5 || log-likelihood : -2270.13551918563
## EM: Iteration : 6 || log-likelihood : -2269.81801187087
## EM: Iteration : 7 || log-likelihood : -2269.21646549831
## EM: Iteration : 8 || log-likelihood : -2268.76855045316
## EM: Iteration : 9 || log-likelihood : -2268.66239434338
## EM: Iteration : 10 || log-likelihood : -2268.5173684659
## EM: Iteration : 11 || log-likelihood : -2267.91195069571
## EM: Iteration : 12 || log-likelihood : -2266.01773290908
## EM: Iteration : 13 || log-likelihood : -2265.54002211479
## EM: Iteration : 14 || log-likelihood : -2265.53997367785
# mhmmr is a ModelMHMMR object. It contains some methods such as 'summary' and 'plot'
mhmmr$summary()
## ----------------------
## Fitted MHMMR model
## ----------------------
## 
## MHMMR model with K = 10 regimes
## 
##  log-likelihood  nu      AIC       BIC
##        -2265.54 349 -2614.54 -3168.778
## 
## Clustering table:
##  1  2  3  4  5  6  7  8  9 10 
## 17 17 17 15 19 19 15 16 18 24 
## 
## 
## ------------------
## Regime 1 (K = 1):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   1.933309e-05 2.124384e-05 1.916561e-05 1.700387e-05 0.0000155451
## X^1 3.821441e-02 4.199128e-02 3.788337e-02 3.361039e-02 0.0307269546
## 
## Covariance matrix:
##                                                  
##   61.88105  64.87395  65.82106  60.31517 104.9246
##   64.87395  76.94748  76.67896  73.18855 121.0847
##   65.82106  76.67896  79.28542  72.55781 123.0497
##   60.31517  73.18855  72.55781  78.18367 114.7617
##  104.92456 121.08471 123.04968 114.76175 198.0223
## ------------------
## Regime 2 (K = 2):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   2.827447e-05 2.754261e-05 2.708866e-05 2.725869e-05 2.777031e-05
## X^1 5.600808e-02 5.455835e-02 5.365913e-02 5.399596e-02 5.500939e-02
## 
## Covariance matrix:
##                                                 
##  480.32550 52.34772 55.92119 100.51901 140.47317
##   52.34772 20.95463 20.26699  30.00160  32.34991
##   55.92119 20.26699 29.17683  37.77143  40.12183
##  100.51901 30.00160 37.77143  61.44997  67.31065
##  140.47317 32.34991 40.12183  67.31065  89.31451
## ------------------
## Regime 3 (K = 3):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   3.558194e-05 3.246503e-05 3.149696e-05 3.230926e-05 3.319651e-05
## X^1 7.063436e-02 6.444691e-02 6.252519e-02 6.413771e-02 6.589899e-02
## 
## Covariance matrix:
##                                                   
##   88.99223 -54.86067 -28.09875 -48.39661 -90.23697
##  -54.86067  55.81858  28.49121  52.40100  76.35080
##  -28.09875  28.49121  16.69622  28.30399  40.00428
##  -48.39661  52.40100  28.30399  55.28000  72.54715
##  -90.23697  76.35080  40.00428  72.54715 112.89830
## ------------------
## Regime 4 (K = 4):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   2.900766e-05 0.0000366562 3.235088e-05 3.385254e-05 5.177285e-05
## X^1 5.770013e-02 0.0729141096 6.435024e-02 6.733725e-02 1.029832e-01
## 
## Covariance matrix:
##                                                       
##  111.026581  1.6044571  3.8107743  2.9773311 -41.89270
##    1.604457 10.2219346  3.9745888  0.5408306  71.31540
##    3.810774  3.9745888  1.8145033  0.4108295  26.25097
##    2.977331  0.5408306  0.4108295  1.1041686  -6.91685
##  -41.892696 71.3154049 26.2509674 -6.9168503 678.07022
## ------------------
## Regime 5 (K = 5):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   3.501514e-05 4.201636e-05 3.672467e-05 3.731638e-05 5.374788e-05
## X^1 6.979823e-02 8.375429e-02 7.320596e-02 7.438546e-02 1.071396e-01
## 
## Covariance matrix:
##                                                   
##   46.60667  60.52686  63.97621  59.44579 -65.48396
##   60.52686  81.17453  84.81540  77.71562 -86.10664
##   63.97621  84.81540  89.75605  82.33438 -92.05009
##   59.44579  77.71562  82.33438  77.44808 -81.83248
##  -65.48396 -86.10664 -92.05009 -81.83248 105.93579
## ------------------
## Regime 6 (K = 6):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   4.115695e-05 5.006855e-05 4.352742e-05 4.632426e-05 5.786691e-05
## X^1 8.223606e-02 1.000424e-01 8.697252e-02 9.256091e-02 1.156244e-01
## 
## Covariance matrix:
##                                                  
##   38.08685  57.74382  41.41296  77.19636 155.3075
##   57.74382  90.62697  64.68069 118.33223 236.1848
##   41.41296  64.68069  46.46492  85.36845 171.5960
##   77.19636 118.33223  85.36845 162.03042 331.5313
##  155.30748 236.18483 171.59596 331.53126 692.8297
## ------------------
## Regime 7 (K = 7):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   4.883203e-05 5.832711e-05 5.059906e-05 5.909917e-05 9.155577e-05
## X^1 9.778129e-02 1.167942e-01 1.013196e-01 1.183402e-01 1.833313e-01
## 
## Covariance matrix:
##                                                  
##  182.3675 113.07680 121.03686  245.4053  748.5601
##  113.0768  75.15479  78.40509  153.5705  469.9886
##  121.0369  78.40509  82.97036  164.3159  502.3604
##  245.4053 153.57046 164.31589  333.5739 1017.5771
##  748.5601 469.98862 502.36041 1017.5771 3107.8847
## ------------------
## Regime 8 (K = 8):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   6.638267e-05 7.103192e-05 6.141522e-05 9.540986e-05 0.0001449247
## X^1 1.331803e-01 1.425078e-01 1.232143e-01 1.914161e-01 0.2907553427
## 
## Covariance matrix:
##                                                     
##  312.7745 295.938209 188.82385  564.0075  184.034515
##  295.9382 292.286889 180.60933  451.1846    1.650569
##  188.8239 180.609325 114.63278  333.2994   94.926965
##  564.0075 451.184554 333.29943 1729.5584 1799.790899
##  184.0345   1.650569  94.92697 1799.7909 3179.904983
## ------------------
## Regime 9 (K = 9):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4) Beta(d = 5)
## 1   6.933623e-05 7.073846e-05 6.096503e-05 6.532265e-05 9.97045e-05
## X^1 1.394006e-01 1.422198e-01 1.225703e-01 1.313313e-01 2.00456e-01
## 
## Covariance matrix:
##                                                      
##   6.859675  -5.77749   1.560809   3.991844   3.031998
##  -5.777490 104.78908  46.212724 287.948580 170.547741
##   1.560809  46.21272  23.598190 142.584729  84.044967
##   3.991844 287.94858 142.584729 907.838460 534.848652
##   3.031998 170.54774  84.044967 534.848652 336.633852
## ------------------
## Regime 10 (K = 10):
## 
## Regression coefficients:
## 
##      Beta(d = 1)  Beta(d = 2)  Beta(d = 3)  Beta(d = 4)  Beta(d = 5)
## 1   7.623422e-05 7.285109e-05 6.514573e-05 0.0000847447 0.0001370561
## X^1 1.536693e-01 1.468497e-01 1.313176e-01 0.1708242987 0.2762711039
## 
## Covariance matrix:
##                                                  
##  139.9193  192.2679  162.2606  502.0126  830.9939
##  192.2679  289.7827  243.1845  739.3775 1193.0375
##  162.2606  243.1845  206.4592  625.3424 1010.7015
##  502.0126  739.3775  625.3424 1920.5276 3124.5251
##  830.9939 1193.0375 1010.7015 3124.5251 5127.9929
mhmmr$plot()

# mhmmr has also two fields, stat and param which are reference classes as well
# Log-likelihood: 

mhmmr$stat$loglik
## [1] -2265.54
# Parameters of the polynomial regressions: 
mhmmr$param$beta
## , , 1
## 
##              [,1]         [,2]         [,3]         [,4]         [,5]
## [1,] 1.933309e-05 2.124384e-05 1.916561e-05 1.700387e-05 0.0000155451
## [2,] 3.821441e-02 4.199128e-02 3.788337e-02 3.361039e-02 0.0307269546
## 
## , , 2
## 
##              [,1]         [,2]         [,3]         [,4]         [,5]
## [1,] 2.827447e-05 2.754261e-05 2.708866e-05 2.725869e-05 2.777031e-05
## [2,] 5.600808e-02 5.455835e-02 5.365913e-02 5.399596e-02 5.500939e-02
## 
## , , 3
## 
##              [,1]         [,2]         [,3]         [,4]         [,5]
## [1,] 3.558194e-05 3.246503e-05 3.149696e-05 3.230926e-05 3.319651e-05
## [2,] 7.063436e-02 6.444691e-02 6.252519e-02 6.413771e-02 6.589899e-02
## 
## , , 4
## 
##              [,1]         [,2]         [,3]         [,4]         [,5]
## [1,] 2.900766e-05 0.0000366562 3.235088e-05 3.385254e-05 5.177285e-05
## [2,] 5.770013e-02 0.0729141096 6.435024e-02 6.733725e-02 1.029832e-01
## 
## , , 5
## 
##              [,1]         [,2]         [,3]         [,4]         [,5]
## [1,] 3.501514e-05 4.201636e-05 3.672467e-05 3.731638e-05 5.374788e-05
## [2,] 6.979823e-02 8.375429e-02 7.320596e-02 7.438546e-02 1.071396e-01
## 
## , , 6
## 
##              [,1]         [,2]         [,3]         [,4]         [,5]
## [1,] 4.115695e-05 5.006855e-05 4.352742e-05 4.632426e-05 5.786691e-05
## [2,] 8.223606e-02 1.000424e-01 8.697252e-02 9.256091e-02 1.156244e-01
## 
## , , 7
## 
##              [,1]         [,2]         [,3]         [,4]         [,5]
## [1,] 4.883203e-05 5.832711e-05 5.059906e-05 5.909917e-05 9.155577e-05
## [2,] 9.778129e-02 1.167942e-01 1.013196e-01 1.183402e-01 1.833313e-01
## 
## , , 8
## 
##              [,1]         [,2]         [,3]         [,4]         [,5]
## [1,] 6.638267e-05 7.103192e-05 6.141522e-05 9.540986e-05 0.0001449247
## [2,] 1.331803e-01 1.425078e-01 1.232143e-01 1.914161e-01 0.2907553427
## 
## , , 9
## 
##              [,1]         [,2]         [,3]         [,4]        [,5]
## [1,] 6.933623e-05 7.073846e-05 6.096503e-05 6.532265e-05 9.97045e-05
## [2,] 1.394006e-01 1.422198e-01 1.225703e-01 1.313313e-01 2.00456e-01
## 
## , , 10
## 
##              [,1]         [,2]         [,3]         [,4]         [,5]
## [1,] 7.623422e-05 7.285109e-05 6.514573e-05 0.0000847447 0.0001370561
## [2,] 1.536693e-01 1.468497e-01 1.313176e-01 0.1708242987 0.2762711039