Skip to content

Commit

Permalink
differences for PR #6
Browse files Browse the repository at this point in the history
  • Loading branch information
@actions-user
actions-user committed Nov 5, 2024
1 parent ca3a6e4 commit 7a1090b
Show file tree
Hide file tree
Showing 9 changed files with 108 additions and 2,615 deletions.
85 changes: 0 additions & 85 deletions config.yaml
View file

This file was deleted.

Binary file modified fig/other-topics-rendered-unnamed-chunk-14-1.png
View file
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file modified fig/other-topics-rendered-unnamed-chunk-15-1.png
View file
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file modified fig/other-topics-rendered-unnamed-chunk-18-1.png
View file
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
20 changes: 10 additions & 10 deletions gaussian-processes.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -196,8 +196,8 @@ gp_samples <- rstan::sampling(gp_model,
SAMPLING FOR MODEL 'anon_model' NOW (CHAIN 1).
Chain 1:
Chain 1: Gradient evaluation took 8.4e-05 seconds
Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.84 seconds.
Chain 1: Gradient evaluation took 8.1e-05 seconds
Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.81 seconds.
Chain 1: Adjust your expectations accordingly!
Chain 1:
Chain 1:
Expand All @@ -214,9 +214,9 @@ Chain 1: Iteration: 800 / 1000 [ 80%] (Sampling)
Chain 1: Iteration: 900 / 1000 [ 90%] (Sampling)
Chain 1: Iteration: 1000 / 1000 [100%] (Sampling)
Chain 1:
Chain 1: Elapsed Time: 32.324 seconds (Warm-up)
Chain 1: 36.752 seconds (Sampling)
Chain 1: 69.076 seconds (Total)
Chain 1: Elapsed Time: 31.544 seconds (Warm-up)
Chain 1: 35.618 seconds (Sampling)
Chain 1: 67.162 seconds (Total)
Chain 1:
```

Expand Down Expand Up @@ -639,8 +639,8 @@ logistic_gp_fit2 <- rstan::sampling(logistic_gp_model2,
SAMPLING FOR MODEL 'anon_model' NOW (CHAIN 1).
Chain 1:
Chain 1: Gradient evaluation took 0.001319 seconds
Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 13.19 seconds.
Chain 1: Gradient evaluation took 0.001301 seconds
Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 13.01 seconds.
Chain 1: Adjust your expectations accordingly!
Chain 1:
Chain 1:
Expand All @@ -657,9 +657,9 @@ Chain 1: Iteration: 1600 / 2000 [ 80%] (Sampling)
Chain 1: Iteration: 1800 / 2000 [ 90%] (Sampling)
Chain 1: Iteration: 2000 / 2000 [100%] (Sampling)
Chain 1:
Chain 1: Elapsed Time: 17.728 seconds (Warm-up)
Chain 1: 17.74 seconds (Sampling)
Chain 1: 35.468 seconds (Total)
Chain 1: Elapsed Time: 17.43 seconds (Warm-up)
Chain 1: 17.379 seconds (Sampling)
Chain 1: 34.809 seconds (Total)
Chain 1:
```
Expand Down
36 changes: 18 additions & 18 deletions md5sum.txt
View file
Original file line number Diff line number Diff line change
@@ -1,19 +1,19 @@
"file" "checksum" "built" "date"
"CODE_OF_CONDUCT.md" "c93c83c630db2fe2462240bf72552548" "site/built/CODE_OF_CONDUCT.md" "2024-09-24"
"LICENSE.md" "b24ebbb41b14ca25cf6b8216dda83e5f" "site/built/LICENSE.md" "2024-09-24"
"config.yaml" "e8f6119a4b795e12f6ac056d6c04edbc" "site/built/config.yaml" "2024-09-24"
"index.md" "a02c9c785ed98ddd84fe3d34ddb12fcd" "site/built/index.md" "2024-09-24"
"links.md" "8184cf4149eafbf03ce8da8ff0778c14" "site/built/links.md" "2024-09-24"
"episodes/bayesian-statistics.Rmd" "03bfec637d0f05fc298b205854267106" "site/built/bayesian-statistics.md" "2024-09-24"
"episodes/stan.Rmd" "46fcb1a4d0ee08863fe6bfa98b74a6f3" "site/built/stan.md" "2024-09-24"
"episodes/mcmc.Rmd" "7ef41fb8c7230bd38e5f558d119cef42" "site/built/mcmc.md" "2024-09-24"
"episodes/hierarchical-models.Rmd" "e8bc0fcc721c2848bb83d52f0fd0c81b" "site/built/hierarchical-models.md" "2024-09-24"
"episodes/model-comparison.Rmd" "541011664e9d0188cfb4b8c4be2629bf" "site/built/model-comparison.md" "2024-09-24"
"episodes/gaussian-processes.Rmd" "98d071cef9e429fe5a6e4b792623707f" "site/built/gaussian-processes.md" "2024-09-24"
"episodes/other-topics.Rmd" "530deb58f400fe9dde55d3d668f821e2" "site/built/other-topics.md" "2024-09-24"
"episodes/exercises.Rmd" "23fc7c68347dcac44b40e3b7743261cf" "site/built/exercises.md" "2024-09-24"
"instructors/instructor-notes.md" "cae72b6712578d74a49fea7513099f8c" "site/built/instructor-notes.md" "2024-09-24"
"learners/reference.md" "1c7cc4e229304d9806a13f69ca1b8ba4" "site/built/reference.md" "2024-09-24"
"learners/setup.md" "bff6e92b34f86fe7610bfb20128e536e" "site/built/setup.md" "2024-09-24"
"profiles/learner-profiles.md" "60b93493cf1da06dfd63255d73854461" "site/built/learner-profiles.md" "2024-09-24"
"renv/profiles/lesson-requirements/renv.lock" "5f816f07f9adc22b61575f81c60df20d" "site/built/renv.lock" "2024-09-24"
"CODE_OF_CONDUCT.md" "c93c83c630db2fe2462240bf72552548" "site/built/CODE_OF_CONDUCT.md" "2024-11-05"
"LICENSE.md" "b24ebbb41b14ca25cf6b8216dda83e5f" "site/built/LICENSE.md" "2024-11-05"
"config.yaml" "e8f6119a4b795e12f6ac056d6c04edbc" "site/built/config.yaml" "2024-11-05"
"index.md" "a02c9c785ed98ddd84fe3d34ddb12fcd" "site/built/index.md" "2024-11-05"
"links.md" "8184cf4149eafbf03ce8da8ff0778c14" "site/built/links.md" "2024-11-05"
"episodes/bayesian-statistics.Rmd" "03bfec637d0f05fc298b205854267106" "site/built/bayesian-statistics.md" "2024-11-05"
"episodes/stan.Rmd" "46fcb1a4d0ee08863fe6bfa98b74a6f3" "site/built/stan.md" "2024-11-05"
"episodes/mcmc.Rmd" "7ef41fb8c7230bd38e5f558d119cef42" "site/built/mcmc.md" "2024-11-05"
"episodes/hierarchical-models.Rmd" "e8bc0fcc721c2848bb83d52f0fd0c81b" "site/built/hierarchical-models.md" "2024-11-05"
"episodes/model-comparison.Rmd" "541011664e9d0188cfb4b8c4be2629bf" "site/built/model-comparison.md" "2024-11-05"
"episodes/gaussian-processes.Rmd" "98d071cef9e429fe5a6e4b792623707f" "site/built/gaussian-processes.md" "2024-11-05"
"episodes/other-topics.Rmd" "530deb58f400fe9dde55d3d668f821e2" "site/built/other-topics.md" "2024-11-05"
"episodes/exercises.Rmd" "23fc7c68347dcac44b40e3b7743261cf" "site/built/exercises.md" "2024-11-05"
"instructors/instructor-notes.md" "cae72b6712578d74a49fea7513099f8c" "site/built/instructor-notes.md" "2024-11-05"
"learners/reference.md" "1c7cc4e229304d9806a13f69ca1b8ba4" "site/built/reference.md" "2024-11-05"
"learners/setup.md" "bff6e92b34f86fe7610bfb20128e536e" "site/built/setup.md" "2024-11-05"
"profiles/learner-profiles.md" "60b93493cf1da06dfd63255d73854461" "site/built/learner-profiles.md" "2024-11-05"
"renv/profiles/lesson-requirements/renv.lock" "822cae49cd98518f55ba293994585da0" "site/built/renv.lock" "2024-11-05"
92 changes: 70 additions & 22 deletions other-topics.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -456,10 +456,10 @@ Formula: time | cens(1 - status) ~ sex + age + ph.karno_status
Regression Coefficients:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept 1.27 0.69 -0.15 2.58 1.00 5012 2490
sex -0.51 0.17 -0.83 -0.19 1.00 4789 2922
age 0.01 0.01 -0.01 0.03 1.00 5207 2779
ph.karno_statushigh -0.36 0.18 -0.72 -0.01 1.00 4534 2669
Intercept 1.27 0.68 -0.11 2.57 1.00 4304 2666
sex -0.50 0.16 -0.83 -0.19 1.00 4768 2842
age 0.01 0.01 -0.01 0.03 1.00 4607 2784
ph.karno_statushigh -0.36 0.18 -0.70 -0.02 1.00 4471 2616
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
Expand All @@ -479,10 +479,10 @@ exp(sum_cox$fixed[,1:4])

``` output
Estimate Est.Error l-95% CI u-95% CI
Intercept 3.5486294 1.985898 0.8641319 13.1790893
sex 0.6031900 1.179877 0.4380762 0.8288425
age 1.0131344 1.009490 0.9948428 1.0323407
ph.karno_statushigh 0.6964744 1.194624 0.4884050 0.9913700
Intercept 3.5576154 1.970680 0.8962140 13.1193479
sex 0.6044998 1.179035 0.4369700 0.8279120
age 1.0130229 1.009435 0.9945742 1.0318059
ph.karno_statushigh 0.6956990 1.193047 0.4959062 0.9847118
```

``` r
Expand Down Expand Up @@ -530,7 +530,7 @@ brms::stancode(fit_cox)
```

``` output
// generated with brms 2.21.0
// generated with brms 2.22.0
functions {
/* distribution functions of the Cox proportional hazards model
* parameterize hazard(t) = baseline(t) * mu
Expand All @@ -545,35 +545,68 @@ functions {
real cox_lhaz(real y, real mu, real bhaz, real cbhaz) {
return log(bhaz) + log(mu);
}
vector cox_lhaz(vector y, vector mu, vector bhaz, vector cbhaz) {
return log(bhaz) + log(mu);
}
// equivalent to the log survival function
real cox_lccdf(real y, real mu, real bhaz, real cbhaz) {
// equivalent to the log survival function
return - cbhaz * mu;
}
real cox_lccdf(vector y, vector mu, vector bhaz, vector cbhaz) {
return - dot_product(cbhaz, mu);
}
real cox_lcdf(real y, real mu, real bhaz, real cbhaz) {
return log1m_exp(cox_lccdf(y | mu, bhaz, cbhaz));
}
real cox_lcdf(vector y, vector mu, vector bhaz, vector cbhaz) {
return sum(log1m_exp(- cbhaz .* mu));
}
real cox_lpdf(real y, real mu, real bhaz, real cbhaz) {
return cox_lhaz(y, mu, bhaz, cbhaz) + cox_lccdf(y | mu, bhaz, cbhaz);
}
real cox_lpdf(vector y, vector mu, vector bhaz, vector cbhaz) {
return sum(cox_lhaz(y, mu, bhaz, cbhaz)) + cox_lccdf(y | mu, bhaz, cbhaz);
}
// Distribution functions of the Cox model in log parameterization
real cox_log_lhaz(real y, real log_mu, real bhaz, real cbhaz) {
return log(bhaz) + log_mu;
}
vector cox_log_lhaz(vector y, vector log_mu, vector bhaz, vector cbhaz) {
return log(bhaz) + log_mu;
}
real cox_log_lccdf(real y, real log_mu, real bhaz, real cbhaz) {
return - cbhaz * exp(log_mu);
}
real cox_log_lccdf(vector y, vector log_mu, vector bhaz, vector cbhaz) {
return - dot_product(cbhaz, exp(log_mu));
}
real cox_log_lcdf(real y, real log_mu, real bhaz, real cbhaz) {
return log1m_exp(cox_log_lccdf(y | log_mu, bhaz, cbhaz));
}
real cox_log_lcdf(vector y, vector log_mu, vector bhaz, vector cbhaz) {
return sum(log1m_exp(- cbhaz .* exp(log_mu)));
}
real cox_log_lpdf(real y, real log_mu, real bhaz, real cbhaz) {
return cox_log_lhaz(y, log_mu, bhaz, cbhaz) +
cox_log_lccdf(y | log_mu, bhaz, cbhaz);
}
real cox_log_lpdf(vector y, vector log_mu, vector bhaz, vector cbhaz) {
return sum(cox_log_lhaz(y, log_mu, bhaz, cbhaz)) +
cox_log_lccdf(y | log_mu, bhaz, cbhaz);
}
}
data {
int<lower=1> N; // total number of observations
vector[N] Y; // response variable
array[N] int<lower=-1,upper=2> cens; // indicates censoring
// censoring indicator: 0 = event, 1 = right, -1 = left, 2 = interval censored
array[N] int<lower=-1,upper=2> cens;
int<lower=1> K; // number of population-level effects
matrix[N, K] X; // population-level design matrix
int<lower=1> Kc; // number of population-level effects after centering
Expand All @@ -588,8 +621,28 @@ data {
int prior_only; // should the likelihood be ignored?
}
transformed data {
// indices of censored data
int Nevent = 0;
int Nrcens = 0;
int Nlcens = 0;
array[N] int Jevent;
array[N] int Jrcens;
array[N] int Jlcens;
matrix[N, Kc] Xc; // centered version of X without an intercept
vector[Kc] means_X; // column means of X before centering
// collect indices of censored data
for (n in 1:N) {
if (cens[n] == 0) {
Nevent += 1;
Jevent[Nevent] = n;
} else if (cens[n] == 1) {
Nrcens += 1;
Jrcens[Nrcens] = n;
} else if (cens[n] == -1) {
Nlcens += 1;
Jlcens[Nlcens] = n;
}
}
for (i in 2:K) {
means_X[i - 1] = mean(X[, i]);
Xc[, i - 1] = X[, i] - means_X[i - 1];
Expand All @@ -598,7 +651,8 @@ transformed data {
parameters {
vector[Kc] b; // regression coefficients
real Intercept; // temporary intercept for centered predictors
simplex[Kbhaz] sbhaz; // baseline coefficients
// baseline hazard coefficients
simplex[Kbhaz] sbhaz;
}
transformed parameters {
real lprior = 0; // prior contributions to the log posterior
Expand All @@ -616,16 +670,10 @@ model {
// initialize linear predictor term
vector[N] mu = rep_vector(0.0, N);
mu += Intercept + Xc * b;
for (n in 1:N) {
// special treatment of censored data
if (cens[n] == 0) {
target += cox_log_lpdf(Y[n] | mu[n], bhaz[n], cbhaz[n]);
} else if (cens[n] == 1) {
target += cox_log_lccdf(Y[n] | mu[n], bhaz[n], cbhaz[n]);
} else if (cens[n] == -1) {
target += cox_log_lcdf(Y[n] | mu[n], bhaz[n], cbhaz[n]);
}
}
// vectorized log-likelihood contributions of censored data
target += cox_log_lpdf(Y[Jevent[1:Nevent]] | mu[Jevent[1:Nevent]], bhaz[Jevent[1:Nevent]], cbhaz[Jevent[1:Nevent]]);
target += cox_log_lccdf(Y[Jrcens[1:Nrcens]] | mu[Jrcens[1:Nrcens]], bhaz[Jrcens[1:Nrcens]], cbhaz[Jrcens[1:Nrcens]]);
target += cox_log_lcdf(Y[Jlcens[1:Nlcens]] | mu[Jlcens[1:Nlcens]], bhaz[Jlcens[1:Nlcens]], cbhaz[Jlcens[1:Nlcens]]);
}
// priors including constants
target += lprior;
Expand Down
Loading

0 comments on commit 7a1090b

Please sign in to comment.