123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243 |
- functions {
- real recs_priors_lpmf(array[] int children,
- int start, int end,
- int n_recs,
- int n_classes,
- real recs_duration,
- array [] real age,
- matrix truth_vocs,
- vector mu_pop_level,
- matrix mu_child_level,
- vector alpha_child_level,
- vector child_dev_age,
- real beta_dev
- ) {
- real ll = 0;
-
- for (k in start:end) {
- real chi_mu = mu_pop_level[1]*exp(
- child_dev_age[children[k-start+1]]*age[k]/12.0/10.0+beta_dev*(mu_child_level[children[k-start+1],2]+mu_child_level[children[k-start+1],3]-mu_pop_level[3]-mu_pop_level[4])*age[k]/12.0/10.0
- );
- ll += gamma_lpdf(
- truth_vocs[k,1]/1000/recs_duration | alpha_child_level[1], alpha_child_level[1]/chi_mu
- );
- ll += gamma_lpdf(
- truth_vocs[k,2:]/1000/recs_duration | alpha_child_level[2:], alpha_child_level[2:]./mu_child_level[children[k-start+1],:]' //'
- );
- }
- return ll;
- }
- }
- // TODO
- // use speech rates to set priors on truth_vocs
- data {
- int<lower=1> n_classes; // number of classes
- // analysis data block
- int<lower=1> n_recs;
- int<lower=1> n_children;
- array[n_recs] int<lower=1> children;
- array[n_recs] real<lower=1> age;
- array[n_recs] int<lower=-1> siblings;
- array[n_recs, n_classes] int<lower=0> vocs;
- array[n_children] int<lower=1> corpus;
- real<lower=0> recs_duration;
- int<lower=1> n_corpora;
- // actual speech rates
- int<lower=1> n_rates;
- int<lower=1> n_speech_rate_children;
- array [n_rates,n_classes] int<lower=0> speech_rates;
- array [n_rates] int group_corpus;
- array [n_rates] real<lower=0> durations;
- array [n_rates] real<lower=0> speech_rate_age;
- array [n_rates] int<lower=-1> speech_rate_siblings;
- array [n_rates] int<lower=1,upper=n_speech_rate_children> speech_rate_child;
- // parallel processing
- int<lower=1> threads;
- }
- transformed data {
- matrix[n_recs, n_classes] vocs_matrix = to_matrix(vocs);
- array[n_speech_rate_children] int<lower=1> speech_rate_child_corpus;
- array[n_children] int<lower=-1> child_siblings;
- array[n_speech_rate_children] int<lower=-1> speech_rate_child_siblings;
- int no_siblings = 0;
- int has_siblings = 0;
- for (k in 1:n_rates) {
- speech_rate_child_corpus[speech_rate_child[k]] = group_corpus[k];
- }
- for (k in 1:n_recs) {
- child_siblings[children[k]] = siblings[k];
- }
- for (c in 1:n_children) {
- if (child_siblings[c] == 0) {
- no_siblings += 1;
- }
- else if (child_siblings[c] > 0) {
- has_siblings += 1;
- }
- }
- for (k in 1:n_rates) {
- speech_rate_child_siblings[speech_rate_child[k]] = speech_rate_siblings[k];
- }
- }
- parameters {
- matrix<lower=0>[n_children,n_classes-1] mu_child_level;
- vector [n_children] child_dev_age;
- // speech rates
- vector<lower=0>[n_classes] alpha_child_level; // variance across recordings for a given child
- array[2] matrix<lower=0>[n_classes-1,n_corpora] alpha_corpus_level; // variance among children
- matrix<lower=0>[n_classes-1,n_corpora] mu_corpus_level; // child-level average
- vector<lower=0>[n_classes-1] alpha_pop_level; // variance among corpora
- vector<lower=0>[n_classes] mu_pop_level; // population level averages
- vector<lower=0>[n_classes-1] alpha_pop;
- matrix<lower=0>[n_classes,n_rates] speech_rate; // truth speech rates observed in annotated clips
- matrix<lower=0>[n_speech_rate_children,n_classes-1] speech_rate_child_level; // expected speech rate at the child-level
- // siblings
- real beta_sib_och; // effect of n of siblings on OCH speech
- real beta_sib_adu; // effect of n of siblings on ADU speech
- real<lower=0,upper=1> p_sib; // prob of having siblings
- vector [n_speech_rate_children] child_dev_speech_age;
- // average effect of age
- real alpha_dev;
- real<lower=0> sigma_dev;
- // effect of excess ADU input
- real beta_dev;
- }
- model {
- // priors on actual speech
- target += reduce_sum(
- recs_priors_lpmf, children, 1,
- n_recs, n_classes, recs_duration, age,
- vocs_matrix,
- mu_pop_level, mu_child_level, alpha_child_level,
- child_dev_age, beta_dev
- );
- vector [2] ll;
- int distrib;
- for (c in 1:n_children) {
- // if there is sibling data
- if (child_siblings[c]>=0) {
- distrib = child_siblings[c]>0?2:1;
- mu_child_level[c,1] ~ gamma(
- alpha_corpus_level[distrib,1,corpus[c]],
- (alpha_corpus_level[distrib,1,corpus[c]]/(mu_corpus_level[1,corpus[c]]*exp(
- child_siblings[c]>0?beta_sib_och:0
- )))
- );
- mu_child_level[c,2:] ~ gamma(
- alpha_corpus_level[distrib,2:,corpus[c]],
- (alpha_corpus_level[distrib,2:,corpus[c]]./mu_corpus_level[2:,corpus[c]]*exp(
- child_siblings[c]>0?beta_sib_adu:0
- ))
- );
- }
- // otherwise
- else {
- // assuming no sibling
- ll[1] = log(p_sib)+gamma_lpdf(
- mu_child_level[c,1] | alpha_corpus_level[2,1,corpus[c]], alpha_corpus_level[2,1,corpus[c]]/(mu_corpus_level[1,corpus[c]]*exp(beta_sib_och))
- );
- ll[1] += gamma_lpdf(
- mu_child_level[c,2] | alpha_corpus_level[2,2,corpus[c]], alpha_corpus_level[2,2,corpus[c]]/(mu_corpus_level[2,corpus[c]]*exp(beta_sib_adu))
- );
- ll[1] += gamma_lpdf(
- mu_child_level[c,3] | alpha_corpus_level[2,3,corpus[c]], alpha_corpus_level[2,3,corpus[c]]/(mu_corpus_level[3,corpus[c]]*exp(beta_sib_adu))
- );
- // assuming sibling
- ll[2] = log(1-p_sib)+gamma_lpdf(
- mu_child_level[c,1] | alpha_corpus_level[1,1,corpus[c]], alpha_corpus_level[1,1,corpus[c]]/(mu_corpus_level[1,corpus[c]])
- );
- ll[2] += gamma_lpdf(
- mu_child_level[c,2] | alpha_corpus_level[1,2,corpus[c]], alpha_corpus_level[1,2,corpus[c]]/(mu_corpus_level[2,corpus[c]])
- );
- ll[2] += gamma_lpdf(
- mu_child_level[c,3] | alpha_corpus_level[1,3,corpus[c]], alpha_corpus_level[1,3,corpus[c]]/(mu_corpus_level[3,corpus[c]])
- );
- target += log_sum_exp(ll);
- }
- }
- alpha_child_level ~ gamma(2,1);
- // speech rates
- mu_pop_level ~ exponential(4); // 250 vocs/hour
- alpha_pop_level ~ gamma(8, 4); // sd = 0.35 x \mu
- alpha_pop ~ gamma(10, 10);
- for (i in 1:n_classes-1) {
- alpha_corpus_level[1,i,:] ~ gamma(4, 4/alpha_pop[i]);
- alpha_corpus_level[2,i,:] ~ gamma(4, 4/alpha_pop[i]);
- mu_corpus_level[i,:] ~ gamma(alpha_pop_level[i],alpha_pop_level[i]/mu_pop_level[i+1]);
- }
- for (g in 1:n_rates) {
- real chi_mu = mu_pop_level[1]*exp(
- child_dev_speech_age[speech_rate_child[g]]*speech_rate_age[g]/12.0/10.0 + beta_dev*(speech_rate_child_level[speech_rate_child[g],2]+speech_rate_child_level[speech_rate_child[g],3]-mu_pop_level[3]-mu_pop_level[4])*speech_rate_age[g]/12.0/10.0
- );
- speech_rate[1,g] ~ gamma(
- alpha_child_level[1],
- alpha_child_level[1]/chi_mu
- );
- speech_rate[2:,g] ~ gamma(
- alpha_child_level[2:],
- (alpha_child_level[2:]./(speech_rate_child_level[speech_rate_child[g],:]')) //'
- );
- speech_rates[g,:] ~ poisson(speech_rate[:,g]*durations[g]*1000);
- }
- for (c in 1:n_speech_rate_children) {
- distrib = child_siblings[c]>0?2:1;
- speech_rate_child_level[c,1] ~ gamma(
- alpha_corpus_level[distrib,1,speech_rate_child_corpus[c]],
- (alpha_corpus_level[distrib,1,speech_rate_child_corpus[c]]/(mu_corpus_level[1,speech_rate_child_corpus[c]]*exp(
- speech_rate_child_siblings[c]>0?beta_sib_och:0
- )))
- );
- speech_rate_child_level[c,2:] ~ gamma(
- alpha_corpus_level[distrib,2:,speech_rate_child_corpus[c]],
- (alpha_corpus_level[distrib,2:,speech_rate_child_corpus[c]]./(mu_corpus_level[2:,speech_rate_child_corpus[c]]*exp(
- speech_rate_child_siblings[c]>0?beta_sib_adu:0
- )))
- );
- }
- child_dev_age ~ normal(alpha_dev, sigma_dev);
- child_dev_speech_age ~ normal(alpha_dev, sigma_dev);
- has_siblings ~ binomial(has_siblings+no_siblings, p_sib);
- p_sib ~ uniform(0, 1);
- beta_sib_och ~ normal(0, 1);
- beta_sib_adu ~ normal(0, 1);
- alpha_dev ~ normal(0, 1);
- sigma_dev ~ exponential(1);
- beta_dev ~ normal(0, 1);
- }
|