1- using StatisticalRethinking , StanSample
2- using Printf, ParetoSmoothedImportanceSampling
1+ using ParetoSmoothedImportanceSampling , StanSample
2+ using CSV, Printf, StatsPlots
33
44ProjDir = @__DIR__
55
@@ -47,8 +47,7 @@ sm1 = SampleModel("roaches1", roaches1_stan; tmpdir)
4747rc1 = stan_sample (sm1; data)
4848
4949if success (rc1)
50- roaches1_df = read_samples (sm1; output_format= :dataframe )
51- precis (roaches1_df[:, [Symbol (" beta.$i " ) for i in 1 : data. K]])
50+ stan_summary (sm1, true )
5251 nt1 = read_samples (sm1)
5352
5453 # Compute LOO and standard error
@@ -129,8 +128,7 @@ sm2 = SampleModel("roaches2", roaches2_stan; tmpdir)
129128rc2 = stan_sample (sm2; data)
130129
131130if success (rc2)
132- roaches2_df = read_samples (sm2; output_format= :dataframe )
133- precis (roaches2_df[:, [Symbol (" beta.$i " ) for i in 1 : data. K]])
131+ read_summary (sm2, true )
134132 nt2 = read_samples (sm2)
135133
136134 # Compute LOO and standard error
@@ -141,14 +139,8 @@ if success(rc2)
141139 @printf (" >> elpd_loo = %.1f, SE(elpd_loo) = %.1f\n " , elpd_loo2, se_elpd_loo2)
142140
143141 # Check the shape parameter k of the generalized Pareto distribution
144- if all (pk2 .< 0.5 )
145- println (" All Pareto k estimates OK (k < 0.5)." )
146- else
147- pk_good = sum (pk2 .<= 0.5 )
148- pk_ok = length (pk2[pk2 .<= 0.7 ]) - pk_good
149- pk_bad = length (pk2[pk2 .<= 1 ]) - pk_good - pk_ok
150- println ((good= pk_good, ok= pk_ok, bad= pk_bad, very_bad= sum (pk2 .> 1 )))
151- end
142+ pk_qualify (pk2) |> display
143+
152144end
153145
154146begin
0 commit comments