Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Quantile.0.25(R) always equals Quantile.0.75(R) #88

Closed
robchallen opened this issue Mar 17, 2020 · 5 comments · Fixed by #89
Closed

Quantile.0.25(R) always equals Quantile.0.75(R) #88

robchallen opened this issue Mar 17, 2020 · 5 comments · Fixed by #89
Assignees
Labels

Comments

@robchallen
Copy link

library(EpiEstim)
data(Flu2009)
T <- nrow(Flu2009$incidence)
t_start <- seq(2, T-6) # starting at 2 as conditional on the past observations
t_end <- t_start + 6 # adding 6 to get 7-day windows as bounds included in window
res_weekly <- EpiEstim::estimate_R(Flu2009$incidence, 
                         method="parametric_si",
                         config = EpiEstim::make_config(list(
                             t_start = t_start,
                             t_end = t_end,
                             mean_si = 2.6, 
                             std_si = 1.5)))
res_weekly$R

results in :

t_start t_end   Mean(R)     Std(R) Quantile.0.025(R) Quantile.0.05(R) Quantile.0.25(R) Median(R) Quantile.0.75(R)
1        2     8 1.7357977 0.40913143        1.02874370       1.12193325        2.4589724 1.7037612        2.4589724
2        3     9 1.7491678 0.36472669        1.10882231       1.19547993        2.3891206 1.7238839        2.3891206

Other quantiles look OK

@jstockwin jstockwin self-assigned this Mar 17, 2020
jstockwin added a commit that referenced this issue Mar 17, 2020
@jstockwin jstockwin added the bug label Mar 17, 2020
@jstockwin
Copy link
Collaborator

Hey @robchallen, thanks for filing the issue.

I agree this looks wrong and indeed can see an obvious problem with the code.

Are you able to confirm if #89 solves the problem?

jstockwin added a commit that referenced this issue Mar 17, 2020
@robchallen
Copy link
Author

I did a:

devtools::install_github("annecori/EpiEstim", ref = "fix-broken-quantiles")

I also tried

devtools::install_github("annecori/EpiEstim", ref = github_pull("89"))

and reran:

T <- nrow(Flu2009$incidence)
t_start <- seq(2, T-6) # starting at 2 as conditional on the past observations
t_end <- t_start + 6 # adding 6 to get 7-day windows as bounds included in window
res_weekly <- EpiEstim::estimate_R(Flu2009$incidence, 
                          method="parametric_si",
                          config = EpiEstim::make_config(list(
                              t_start = t_start,
                              t_end = t_end,
                              mean_si = 2.6, 
                              std_si = 1.5)))
glimpse(res_weekly$R)

both time result was as follows & problem persists:

Observations: 25
Variables: 11
$ t_start             <dbl> 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26
$ t_end               <dbl> 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32
$ `Mean(R)`           <dbl> 1.7357977, 1.7491678, 1.5370581, 1.4318389, 1.4227247, 1.6353733, 1.6639845, 1.6161467, 1.4344239, 1.4…
$ `Std(R)`            <dbl> 0.40913143, 0.36472669, 0.30741163, 0.27059213, 0.25150457, 0.25234358, 0.23300443, 0.20692639, 0.1765…
$ `Quantile.0.025(R)` <dbl> 1.02874370, 1.10882231, 0.99470299, 0.95144658, 0.97314263, 1.17863324, 1.23894574, 1.23622420, 1.1093…
$ `Quantile.0.05(R)`  <dbl> 1.12193325, 1.19547993, 1.06869352, 1.01766106, 1.03580815, 1.24358968, 1.30015068, 1.29149686, 1.1568…
$ `Quantile.0.25(R)`  <dbl> 2.4589724, 2.3891206, 2.0751763, 1.9040473, 1.8601072, 2.0713724, 2.0648767, 1.9708944, 1.7366694, 1.6…
$ `Median(R)`         <dbl> 1.7037612, 1.7238839, 1.5166133, 1.4148298, 1.4079324, 1.6224126, 1.6531215, 1.6073240, 1.4271858, 1.4…
$ `Quantile.0.75(R)`  <dbl> 2.4589724, 2.3891206, 2.0751763, 1.9040473, 1.8601072, 2.0713724, 2.0648767, 1.9708944, 1.7366694, 1.6…
$ `Quantile.0.95(R)`  <dbl> 2.4589724, 2.3891206, 2.0751763, 1.9040473, 1.8601072, 2.0713724, 2.0648767, 1.9708944, 1.7366694, 1.6…
$ `Quantile.0.975(R)` <dbl> 2.6247813, 2.5331192, 2.1955399, 2.0088487, 1.9563365, 2.1657455, 2.1507413, 2.0461979, 1.8005896, 1.7…

also I notice quantile 0,95 is the same as 0,75 and 0,25. Not sure if that was the case before or not.

@jstockwin
Copy link
Collaborator

Hm, my fix should solve those other quantiles too - I did notice they were also wrong.

It could just be some install thing your side, or perhaps my fix doesn't work. I'll try and find time to take a proper look at some point...

@robchallen
Copy link
Author

There must have been some caching issue when I retested it. Looking at this now the quantile results are appropriate and all different. I think this is fixed.

@jstockwin
Copy link
Collaborator

This fix has been pushed to CRAN, in version 2.2-3.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

Successfully merging a pull request may close this issue.

2 participants