# Finding biological condition-specific changes in T- and B-cell receptor repertoires with ClustIRR

    library(knitr)
    library(ClustIRR)
    library(ggplot2)
    library(patchwork)
    library(ggforce)
    library(rBLAST)
    theme_set(new = theme_bw(base_size = 10))

Note: BLAST was not installed when this vignette was built!

    # System requirements
    #  The BLAST+ software needs to be installed on your system. Installation
    # instructions are available in this package's
    # [INSTALL](https://github.com/mhahsler/rBLAST/blob/devel/INSTALL) file and
    # at \url{https://www.ncbi.nlm.nih.gov/books/NBK569861/}.

    # R needs to be able to find the executable. After installing the software,
    # try in R
    Sys.which("blastp")

    # If the command returns "" instead of the path to the executable,
    # then you need to set the environment variable called PATH. In R
    Sys.setenv(PATH = paste(Sys.getenv("PATH"),
       "path_to_your_BLAST_installation", sep=.Platform$path.sep))

    data("D2", package = "ClustIRR")

    cl <- clustirr(s = D2, control = list(blast_gmi = 0.8))

    gcd <- detect_communities(graph = cl$graph,
                              algorithm = "leiden",
                              metric = "average",
                              resolution = 1,
                              iterations = 100,
                              chains = c("CDR3a", "CDR3b"))

    dim(gcd$community_occupancy_matrix)

    honeycomb <- get_honeycombs(com = gcd$community_occupancy_matrix)
    wrap_plots(honeycomb, nrow = 5, ncol = 3)+
        plot_annotation(tag_levels = 'A')

    d <- dco(community_occupancy_matrix = gcd$community_occupancy_matrix,
             groups = c(1, 1, 1, 2, 2, 2),
             mcmc_control = list(mcmc_warmup = 300,
                                 mcmc_iter = 600,
                                 mcmc_chains = 2,
                                 mcmc_cores = 1,
                                 mcmc_algorithm = "NUTS",
                                 adapt_delta = 0.9,
                                 max_treedepth = 10))

    ggplot(data = d$posterior_summary$beta)+
        geom_sina(aes(x = sample, y = mean))|
        ggplot(data = d$posterior_summary$beta_mu)+
        geom_sina(aes(x = as.character(g), y = mean))

    ggplot(data = d$posterior_summary$epsilon)+
        facet_wrap(facets = ~contrast, ncol = 2)+
        geom_errorbar(aes(x = community, y = mean, ymin = L95, ymax = H95),
                      col = "lightgray", width = 0)+
        geom_point(aes(x = community, y = mean), shape = 21, fill = "black",
                   stroke = 0.4, col = "white", size = 1.25)+
        theme(legend.position = "top")+
        ylab(label = expression(delta))+
        scale_x_continuous(expand = c(0,0))