From 043f02dcd15804d6a51a7ed4273ea2016320b507 Mon Sep 17 00:00:00 2001 From: michellepistner Date: Fri, 12 Jan 2024 13:05:57 -0500 Subject: [PATCH] another small bug --- R/fit_basset.R | 11 ++++---- vignettes/non-linear-models.Rmd | 48 ++++++++++++++++----------------- 2 files changed, 30 insertions(+), 29 deletions(-) diff --git a/R/fit_basset.R b/R/fit_basset.R index b1c1e0c..1b06b2e 100644 --- a/R/fit_basset.R +++ b/R/fit_basset.R @@ -92,11 +92,6 @@ basset <- function(Y=NULL, X, upsilon=NULL, Theta=NULL, Gamma=NULL, Xi=NULL, lin stop("Theta and Gamma must be of the same length.") } - ## setting newdata <- X if newdata is null - if(is.null(newdata)){ - newdata <- X - } - ## evaluating theta and gamma ## theta theta_eval <- function(Theta, X, linear){ @@ -139,6 +134,7 @@ basset <- function(Y=NULL, X, upsilon=NULL, Theta=NULL, Gamma=NULL, Xi=NULL, lin Gamma_comb <- Reduce('+', Gamma_trans) ## fitting the joint model + ## newdata auto handled by pibble collapse_samps <- pibble(Y, X=diag(ncol(X)), upsilon, Theta_comb, Gamma_comb, Xi, init, pars, newdata = newdata, ...) ## fitting uncollapse using the joint samples @@ -176,6 +172,11 @@ basset <- function(Y=NULL, X, upsilon=NULL, Theta=NULL, Gamma=NULL, Xi=NULL, lin } } + ## setting newdata <- X if newdata is null + if(is.null(newdata)){ + newdata <- X + } + for(i in 1:num.comp){ ## if num.comp == 1 --> return the samples from Lambda above if(num.comp == 1){ diff --git a/vignettes/non-linear-models.Rmd b/vignettes/non-linear-models.Rmd index 669a0d3..214c0b3 100644 --- a/vignettes/non-linear-models.Rmd +++ b/vignettes/non-linear-models.Rmd @@ -145,30 +145,30 @@ predicted <- predict(fit.clr, X_predict, jitter=1) Now I am going to create a visual that shows the observed data in CLR coordinates (to do that I will need to add a pseudo-count) along with the smoothed estimates. ```{r fig.height=5, fig.width=7} -# family_names <- as(mallard_family$tax_table$Family, "vector") -# Y_clr_tidy <- clr_array(Y+0.65, parts = 1) %>% -# gather_array(mean, coord, sample) %>% -# mutate(time = X[1,sample], -# coord = paste0("CLR(", family_names[coord],")")) -# -# predicted_tidy <- gather_array(predicted, val, coord, sample, iter) %>% -# mutate(time = X_predict[1,sample]) %>% -# filter(!is.na(val)) %>% -# group_by(time, coord) %>% -# summarise_posterior(val, na.rm=TRUE) %>% -# ungroup() %>% -# mutate(coord = paste0("CLR(", family_names[coord],")")) -# -# ggplot(predicted_tidy, aes(x = time, y=mean)) + -# geom_ribbon(aes(ymin=p2.5, ymax=p97.5), fill="darkgrey", alpha=0.5) + -# geom_ribbon(aes(ymin=p25, ymax=p75), fill="darkgrey", alpha=0.9)+ -# geom_line(color="blue") + -# geom_point(data = Y_clr_tidy, alpha=0.5) + -# facet_wrap(~coord, scales="free_y") + -# theme_minimal()+ -# theme(axis.title.y = element_blank(), -# axis.title.x = element_blank(), -# axis.text.x = element_text(angle=45)) +family_names <- as(mallard_family$tax_table$Family, "vector") +Y_clr_tidy <- clr_array(Y+0.65, parts = 1) %>% + gather_array(mean, coord, sample) %>% + mutate(time = X[1,sample], + coord = paste0("CLR(", family_names[coord],")")) + +predicted_tidy <- gather_array(predicted, val, coord, sample, iter) %>% + mutate(time = X_predict[1,sample]) %>% + filter(!is.na(val)) %>% + group_by(time, coord) %>% + summarise_posterior(val, na.rm=TRUE) %>% + ungroup() %>% + mutate(coord = paste0("CLR(", family_names[coord],")")) + +ggplot(predicted_tidy, aes(x = time, y=mean)) + + geom_ribbon(aes(ymin=p2.5, ymax=p97.5), fill="darkgrey", alpha=0.5) + + geom_ribbon(aes(ymin=p25, ymax=p75), fill="darkgrey", alpha=0.9)+ + geom_line(color="blue") + + geom_point(data = Y_clr_tidy, alpha=0.5) + + facet_wrap(~coord, scales="free_y") + + theme_minimal()+ + theme(axis.title.y = element_blank(), + axis.title.x = element_blank(), + axis.text.x = element_text(angle=45)) ```