Notes for Jacques images.

1 Changelog

• 20230317: Initial cellpose implementation.
• 20230315: Separation of input image into timepoints.
• 20230310: Used geopandas to trace cells over time by position.
• 20230307: Setting up my environment to handle ipynb and python markdown.

2 Settings which may be necessary in emacs

I have been messing with my emacs setup to make this project easier.

(defvar elpy-rpc-virtualenv-path 'current)
(pyvenv-activate "/sw/local/fiji/202211")
(pyvenv-activate "venv")

3 Introduction

Jacques is seeking to segment and follow cells over time.

Steps performed so far:

1. Python/venv installation of fiji
2. Figured out some easy interactions between python and actual image data.
3. Implemented a simple function to find minimum distances between putative cells using the extant methods.

4 Next steps

1. Note that the czi images are immediately openable via fiji’s bioformat interface.
2. Load dataset via fiji
3. Invoke cellpose via the python interface (same as pyimagej)
4. Save ROIs produced by cellpose, save them to zip output file

1. The roi interface in fiji can address these

5. Currently a macro is invoked which performs a set of measurements on every ROI and saves them to a csv. This creates the primary data structure used for all following processes.
6. Need to create a datastructure which identifies each individual cell over time given these ROIs which x/y/time(frame) along with the measurements of interest (area/mean/stdev/intensities by color).
   1. Small postprocessing details: the intensity values produced must be normalized by cell area.

5 Demonstrate the first couple of steps using an actual dataset

Jacques sent me an image acquired from the microscope and I saved it as ‘test_data/raw.tif’, he also sent me a cellpose model which I saved to the ‘models/’ directory.

Given these as a starting point, let us try to open the image with a fiji instance and split it into a series of smaller files by timepoint.

5.1 Load necessary python modules

import os
import glob
import imagej
import pandas
from pandas import DataFrame
import numpy as np
import geopandas
import matplotlib.pyplot as plt
import seaborn
import scyjava
import multiprocessing as mp
from cellpose import models, io
from cellpose.io import *
from collections import defaultdict
from jpype import JArray, JInt

5.2 Set up some initial variables

base_dir = '/lab/scratch/atb/imaging/mtb_2023'
os.chdir(base_dir)
input_file = f"{base_dir}/test_data/raw.tif"
output_dir = f"{base_dir}/outputs"

5.3 Start imagej/fiji

Note that I am using fiji/imagej from within a virtual environment which I symbolically linked to the local directory ‘venv/’.

As a result, when I initialize fiji, I will call the base directory of the downloaded fiji tree within the virtual environment, which I somewhat erroneously put in bin/.

scyjava.config.add_option('-Xmx64g')
start_dir = os.getcwd()
ij = imagej.init('venv/bin/Fiji.app', mode = 'interactive')
## Something about this init() function changes the current working directory.
os.chdir(start_dir)
ij.getVersion()

## '2.9.0/1.53t'

5.4 Open the input file and split it by time

Open the file, figure out its dimensions, and write portions of it to an output directory.

I wrote this thinking I could parallelize the output writing to 8 cpus. But I think I do not yet understand how python scopes variables in this context and so it did not quite work. I turned that off for the moment and ran it and it finished in about a minute.

The following function may require a test to see if the output directory already exists because scijava will freak out.

def separate_slices(input_file, output_base = 'outputs', wanted_x = True, wanted_y = True,
                    wanted_z = 1, wanted_channel = 2, cpus = 8, overwrite = False):
    pool = mp.Pool(cpus)
    input_base = os.path.basename(input_file)
    input_dir = os.path.dirname(input_file)
    input_name = os.path.splitext(input_base)[0]
    output_directory = f"{input_dir}/outputs/{input_name}_substack"
    os.makedirs(output_directory, exist_ok = True)
    print("Starting to open the input file, this takes a moment.")
    dataset = ij.io().open(input_file)
    print(f"Opened input file, writing images to {output_directory}")

    data_info = {}
    for element in range(len(dataset.dims)):
        name = dataset.dims[element]
        data_info[name] = dataset.shape[element]
    print(f"This dataset has dimensions: X:{data_info['X']} Y:{data_info['Y']} Z:{data_info['Z']} Time:{data_info['Time']}")

    def save_slice(timepoint):
        wanted_slice = dataset[:, :, wanted_channel, wanted_z, timepoint]
        slice_image = ij.py.to_dataset(wanted_slice)
        output_filename = f"{output_directory}/{input_name}_{timepoint}.tif"
        if (os.path.exists(output_filename)):
            if overwrite:
                print(f"Rewriting {output_filename}")
                os.remove(output_filename)
                saved = ij.io().save(slice_image, output_filename)
            else:
                print(f"Skipping {output_filename}, it already exists.")
        else:
            saved = ij.io().save(slice_image, output_filename)
            print(f"Saving image {input_name}_{timepoint}.")
        return wanted_slice

    slices = []
    for timepoint in range(data_info['Time']):
        saved = save_slice(timepoint)
        slices.append(saved)
        ## I want to see if I can split the saves across 8 cpus.
        ## pool.apply(save_slice, args = (timepoint, ))
    ## pool.close()
    return dataset, slices, output_directory

raw_split, saved_slices, output_directory = separate_slices(input_file)

5.5 Cellpose

At this point we should have a directory containing files of individual timepoints. Jacques sent me an initial implementation of the usage of cellpose to call individual cells. Let us include that now. I think the previous function should probably also return the directory of the separated input files.

def invoke_cellpose(input_directory, model_file, channels = [[0, 0]], diameter = 160,
                    threshold = 0.4, do_3D = False, batch_size = 64):
    """ Invoke cellpose using the split files from the previous function.

    Jacques mentioned that this is getting slightly different results than the
    GUI-invoked version of this same function.  I figured that since I have a minute,
    I can poke at this and see if I can understand the requisites before he arrives.
    """

    ## Relevant options:
    ## model_type(cyto, nuclei, cyto2), net_avg(T/F if load built in networks and average them)
    model = models.CellposeModel(pretrained_model = model_file)
    files = get_image_files(input_directory, '_masks', look_one_level_down = False)
    imgs = [imread(f) for f in files]
    nimg = len(imgs)
    print(f"Read {nimg} images, starting cellpose.")
    ## Relevant options:
    ## batch_size(increase for more parallelization), channels(two element list of two element
    ## channels to segment; the first is the segment, second is optional nucleus;
    ## internal elements are color channels to query, so [[0,0],[2,3]] means do main cells in
    ## grayscale and a second with cells in blue, nuclei in green.
    ## channel_axis, z_axis ? invert (T/F flip pixels from b/w I assume),
    ## normalize(T/F percentile normalize the data), diameter, do_3d,
    ## anisotropy (rescaling factor for 3d segmentation), net_avg (average models),
    ## augment ?, tile ?, resample, interp, flow_threshold, cellprob_threshold (interesting),
    ## min_size (turned off with -1), stitch_threshold ?, rescale ?.
    output_files = defaultdict(dict)
    existing_files = 0
    for filename in files:
        f_name = os.path.splitext(filename)[0]
        output_mask = f"{start_dir}/{f_name}_cp_masks.png"
        output_txt = f"{start_dir}/{f_name}_cp_outlines.txt"
        output_files[filename]['input_file'] = f"{start_dir}/{filename}"
        output_files[filename]['output_mask'] = output_mask
        output_files[filename]['output_txt'] = output_txt
        output_files[filename]['exists'] = False
        if (os.path.exists(output_files[filename]['output_txt'])):
            existing_files = existing_files + 1
            output_files[filename]['exists'] = True

    if existing_files > 0:
        print(f"Out of {nimg} output files, {existing_files} already exist, not running cellpose.")
    else:
        masks, flows, styles = model.eval(
            imgs, diameter = diameter, channels = channels, flow_threshold = threshold,
            do_3D = do_3D, batch_size = batch_size)
        io.save_to_png(imgs, masks, flows, files)
    print("Finished cellpose, returning output filenames.")
    return output_files

os.chdir(start_dir)
output_directory = 'test_data/outputs/raw_substack'
cellpose_result = invoke_cellpose(output_directory, 'models/CP_20220523_104016')

6 Create Regions of interest from cellpose outputs

In Jacques notebook, it looks like he only extracts ROIs from one of the cellpose slices. I am assuming the goal is to extend this across all images?

There is an important caveat that I missed: imagej comes with a python2-based scripting language from which it appears some of his code is coming. As a result I should look carefully before using it, and pay close attention to the examples provided here for the most appropriate ways of interacting with the ROI manager etc:

https://github.com/imagej/pyimagej/blob/main/doc/examples/blob_detection_interactive.py

# get ImageJ resources
PolygonRoi = scyjava.jimport('ij.gui.PolygonRoi')
Overlay = scyjava.jimport('ij.gui.Overlay')
Regions = scyjava.jimport('net.imglib2.roi.Regions')
LabelRegions = scyjava.jimport('net.imglib2.roi.labeling.LabelRegions')
ov = Overlay()
from pandas import DataFrame
os.chdir(start_dir)

## The following is from a mix of a couple of implementations I found:
## https://pyimagej.readthedocs.io/en/latest/Classic-Segmentation.html
## an alternative method may be taken from:
## https://pyimagej.readthedocs.io/en/latest/Classic-Segmentation.html#segmentation-workflow-with-imagej2
## My goal is to pass the ROI regions to this function and create a similar df.

def slices_to_roi_measurements(cellpose_result, saved_slices):
    output_dict = cellpose_result
    cellpose_slices = list(cellpose_result.keys())
    slice_number = 0
    for slice_name in cellpose_slices:
        output_dict[slice_name]['slice_number'] = slice_number
        slice_data = saved_slices[slice_number]
        input_tif = cellpose_result[slice_name]['input_file']
        input_txt = cellpose_result[slice_name]['output_txt']
        input_mask = cellpose_result[slice_name]['output_mask']
        print(f"Processing cellpose outline: {input_txt}")
        # convert Dataset to ImagePlus
        imp = ij.py.to_imageplus(slice_data)
        rm = ij.RoiManager.getRoiManager()
        ## The logic for this was taken from:
        ## https://stackoverflow.com/questions/73849418/is-there-any-way-to-switch-imagej-macro-code-to-python3-code
        txt_fh = open(input_txt, 'r')
        ij.IJ.run("Set Measurements...", "area mean min centroid median skewness kurtosis redirect=None decimal=3")
        roi_stats = defaultdict(list)
        for line in txt_fh:
            xy = line.rstrip().split(",")
            xy_coords = [int(element) for element in xy]
            x_coords = [int(element) for element in xy[::2]]
            y_coords = [int(element) for element in xy[1::2]]
            xcoords_jint = JArray(JInt)(x_coords)
            ycoords_jint = JArray(JInt)(y_coords)
            polygon_roi_instance = scyjava.jimport('ij.gui.PolygonRoi')
            roi_instance = scyjava.jimport('ij.gui.Roi')
            imported_polygon = polygon_roi_instance(xcoords_jint, ycoords_jint,
                                                    len(x_coords), int(roi_instance.POLYGON))
            imp.setRoi(imported_polygon)
            ij.IJ.run(imp, 'Measure', '')
        slice_result = ij.ResultsTable.getResultsTable()
        slice_table = ij.convert().convert(slice_result, scyjava.jimport('org.scijava.table.Table'))
        slice_measurements = ij.py.from_java(slice_table)
        output_dict[slice_name]['measurements'] = slice_measurements
        ij.IJ.run('Clear Results')
        txt_fh.close()
        imp.setOverlay(ov)
        imp.getProcessor().resetMinAndMax()
        slice_number = slice_number + 1
    return output_dict

slice_measurements = slices_to_roi_measurements(cellpose_result, saved_slices)

7 Convert the slice measurements to pandas df

slices_to_roi_measurements() returns a dictionary with keys which are the filenames of each raw tif file. Each element of that dictionary is in turn a dictionary containing some information about the files along with a df of the measurements provided by imagej.

My little geopandas function assumes a single long df with some columns which tell it which timepoint. So lets make a quick function to give that here. OTOH it may be wiser/better to make some changes to slices_to_roi_measurements() so that it returns that format df; but since I am using this as a learning experience to get more comfortable with python data structures, I will not do it that way.

## def convert_slices_to_pandas(slices):

concatenated = pandas.DataFrame()
slice_keys = list(slice_measurements.keys())
slice_counter = 0
for k in slice_keys:
    slice_counter = slice_counter + 1
    current_slice = slice_measurements[k]
    print(f"The slice is {k}")
    slice_number = current_slice['slice_number']
    slice_data = current_slice['measurements']
    slice_data['slice_number'] = slice_number
    if (slice_counter == 1):
        concatenated = slice_data
    else:
        concatenated = pandas.concat([concatenated, slice_data])

8 Create cell groups

8.1 Load the TSV data

In the above R implementation, I dumped the data to a tsv file, and will just read that here for simplicity.

While I am testing, I will extract only the columns containing the ID, x, y, and time(Frame).

9 Load the coordinates into a geospatial distance method

geopandas looks like it handles coordinate data in a similar fashion to R’s spatstat.geom. I am not sure if it will be able to handle a z-stack; but I assume it will prove easy to extend the implementation to points_from_xyz().

gdf = geopandas.GeoDataFrame(
    concatenated,
    geometry = geopandas.points_from_xy(concatenated.X, concatenated.Y))

I think the easiest way to extract minimum distances in python is to split this df into pairs by Frame followed by using the sjoin_nearest() function to find nearest points from a->b.

9.1 subset the gdf data

In the following block, I will iterate over every time from 1 to end-1 and extract the elements of the gdf from those two time points. Then I will use the sjoin_nearest() to get me every time-pairwise closest distance as a new dataframe.

These pairwise distance dataframes will be appended to the eponymous pairwise_distances array; which upon completion will comprise n-1 elements where the first element includes the results for timepoints 1,2 and the last element includes n-1,n.

final_time = gdf.Frame.max()
pairwise_distances = []
for start_time in range(1, final_time):
    i = start_time
    j = i + 1
    ti_idx = gdf.Frame == i
    tj_idx = gdf.Frame == j
    print("Getting distances of dfs ", i, " and ", j, ".")
    ti = gdf[ti_idx]
    tj = gdf[tj_idx]
    ti_rows = ti.shape[0]
    tj_rows = tj.shape[0]
    titj = geopandas.sjoin_nearest(ti, tj, distance_col = "pairwise_dist")
    pairwise_distances.append(titj)

9.2 Filter the pairwise distances and find neighbors over time

Given the pairwise_distances array, I should be able to iterate over it and make a series of arrays containing the cell IDs which are composed of cells which are closest to each other over time.

The primary caveat here is that we know that sometimes the cells are merged or other shenanigans happen which cause the closest distance to get relatively large, thus I will filter out any pairwise distances which are greater than max_dist.

As a result, the number of cells observed over time will by definition not be the ~ 20 observed in the 121 timepoints, but will be some larger number comprised of the sum over time. A likely next step will be to merge these partial groups post-facto, but not right now, I just want to get this initial implementation working.

Here is the logic:

1. Use id_counter to make a global cell ID across all timepoints.
2. Use traced{} to associate the id_counter cell ID to an array of observed cells over time.
3. Use ends{} to associate every known current endpoint with a id_counter cell ID so that we can keep track of who is who.
4. Then iterate over every element of pairwise_distances and do the following:

1. drop elements with distances > max_dist.
2. check to see if the early-timepoint nearest neighbor (start_cell) has been seen in ends{}: * If so add end_cell to it and to the appropriate element of traced{}. * If not, initialize traced{} with the start/end and add it to ends{}.

Upon completion I should have a dictionary of arrays containing all of the contiugous nearest-neighbor cells that stayed within max_dist of each other.

## Start over if the distance between two cells is bigger than 10 units.
max_dist = 10.0
## Start over if the proportion between the bigger and smaller cell area has changed
## by more than a factor of 0.25.
max_prop = 0.7
id_counter = 0
## Cell IDs pointing to a list of cells
traced = {}
## Endpoints pointing to the cell IDs
ends = {}
for i in range(0, final_time - 1):
    query = pairwise_distances[i]

    passed_idx = query.pairwise_dist <= max_dist
    failed_idx = query.pairwise_dist > max_dist
    if (failed_idx.sum() > 0):
        print("Skipped ", failed_idx.sum(), " elements in segment ", i,
              " because distance is too large.")
        query = query[passed_idx]

    prop_change = query.Area_left / query.Area_right
    increased_idx = prop_change > 1.0
    prop_change[increased_idx] = 1.0 / prop_change[increased_idx]
    failed_idx = prop_change < max_prop
    passed_idx = prop_change >= max_prop
    if (failed_idx.sum() > 0):
        print("Skipped ", failed_idx.sum(), " elements in segment ", i,
              " because the size changed too much.")
        query = query[passed_idx]

    for row in query.itertuples():
        start_cell = row.ID_left
        end_cell = row.ID_right
        if start_cell in ends.keys():
            cell_id = ends[start_cell]
            current_value = traced[cell_id]
            current_value.append(end_cell)
            traced[cell_id] = current_value
            ends[end_cell] = cell_id
        else:
            id_counter = id_counter + 1
            traced[id_counter] = [start_cell, end_cell]
            ends[end_cell] = id_counter

traced

Ok, I looked at the result of traced and it seems to make sense to me; lets rewrite this as a single function that Jacques can add to his notebook. I will assume he has a pandas df as input.

9.3 Get information from a group of cells

As a final step, we should be able to extract and play with the information from one or more groups of cells.

cell_id = 6
cell_idx = input['ID'].isin(traced[cell_id])
cell_data = input[cell_idx].reset_index()

plt.scatter(cell_data['X'], cell_data['Y'])
final_row = cell_data.Frame.max()
for start_time in range(1, final_row - 1):
    ti_idx = cell_data.Frame == start_time
    tj_idx = cell_data.Frame == start_time + 1
    p1x = cell_data[ti_idx].X
    p2x = cell_data[tj_idx].X
    p1y = cell_data[ti_idx].Y
    p2y = cell_data[tj_idx].Y
    x_points = [p1x, p2x]
    y_points = [p1y, p2y]
    plt.plot(x_points, y_points)
plt.show()

seaborn.violinplot(data = cell_data.Area)
plt.show()

10 Single function definition

Now let us combine the things learned above into a single function.

def nearest_cells_over_time(df, max_dist = 10.0, x_column = 'X', y_column = 'Y'):
    """Trace cells over time"""
    gdf = geopandas.GeoDataFrame(
        df,
        geometry = geopandas.points_from_xy(df[x_column], df[y_column]))

    final_time = gdf.Frame.max()
    pairwise_distances = []
    for start_time in range(1, final_time):
        i = start_time
        j = i + 1
        ti_idx = gdf.Frame == i
        tj_idx = gdf.Frame == j
        print("Getting distances of dfs ", i, " and ", j, ".")
        ti = gdf[ti_idx]
        tj = gdf[tj_idx]
        ti_rows = ti.shape[0]
        tj_rows = tj.shape[0]
        titj = geopandas.sjoin_nearest(ti, tj, distance_col = "pairwise_dist")
        pairwise_distances.append(titj)

    max_dist = 10.0
    id_counter = 0
    ## Cell IDs pointing to a list of cells
    traced = {}
    ## Endpoints pointing to the cell IDs
    ends = {}
    for i in range(0, final_time - 1):
        query = pairwise_distances[i]
        passed_idx = query.pairwise_dist <= max_dist
        failed_idx = query.pairwise_dist > max_dist
        if (failed_idx.sum() > 0):
            print("Skipped ", failed_idx.sum(), " elements in segment ", i, ".")
        query = query[passed_idx]
        for row in query.itertuples():
            start_cell = row.ID_left
            end_cell = row.ID_right
            if start_cell in ends.keys():
                cell_id = ends[start_cell]
                current_value = traced[cell_id]
                current_value.append(end_cell)
                traced[cell_id] = current_value
                ends[end_cell] = cell_id
            else:
                id_counter = id_counter + 1
                traced[id_counter] = [start_cell, end_cell]
                ends[end_cell] = id_counter

    return traced

pander::pander(sessionInfo())

R version 4.2.0 (2022-04-22)

Platform: x86_64-pc-linux-gnu (64-bit)

locale: LC_CTYPE=en_US.UTF-8, LC_NUMERIC=C, LC_TIME=en_US.UTF-8, LC_COLLATE=en_US.UTF-8, LC_MONETARY=en_US.UTF-8, LC_MESSAGES=en_US.UTF-8, LC_PAPER=en_US.UTF-8, LC_NAME=en_US.UTF-8, LC_ADDRESS=en_US.UTF-8, LC_TELEPHONE=en_US.UTF-8, LC_MEASUREMENT=en_US.UTF-8 and LC_IDENTIFICATION=en_US.UTF-8

attached base packages: stats4, stats, graphics, grDevices, utils, datasets, methods and base

other attached packages: spatstat.geom(v.3.0-6), spatstat.data(v.3.0-0), hpgltools(v.1.0), testthat(v.3.1.6), reticulate(v.1.28), SummarizedExperiment(v.1.28.0), GenomicRanges(v.1.50.2), GenomeInfoDb(v.1.34.9), IRanges(v.2.32.0), S4Vectors(v.0.36.1), MatrixGenerics(v.1.10.0), matrixStats(v.0.63.0), Biobase(v.2.58.0) and BiocGenerics(v.0.44.0)

loaded via a namespace (and not attached): rappdirs(v.0.3.3), rtracklayer(v.1.58.0), tidyr(v.1.3.0), ggplot2(v.3.4.1), clusterGeneration(v.1.3.7), bit64(v.4.0.5), knitr(v.1.42), DelayedArray(v.0.24.0), data.table(v.1.14.6), KEGGREST(v.1.38.0), RCurl(v.1.98-1.10), doParallel(v.1.0.17), generics(v.0.1.3), GenomicFeatures(v.1.50.4), callr(v.3.7.3), RhpcBLASctl(v.0.23-42), cowplot(v.1.1.1), usethis(v.2.1.6), RSQLite(v.2.2.20), shadowtext(v.0.1.2), bit(v.4.0.5), enrichplot(v.1.18.3), xml2(v.1.3.3), httpuv(v.1.6.8), assertthat(v.0.2.1), viridis(v.0.6.2), xfun(v.0.37), hms(v.1.1.2), jquerylib(v.0.1.4), evaluate(v.0.20), promises(v.1.2.0.1), fansi(v.1.0.4), restfulr(v.0.0.15), progress(v.1.2.2), caTools(v.1.18.2), dbplyr(v.2.3.0), igraph(v.1.4.0), DBI(v.1.1.3), htmlwidgets(v.1.6.1), purrr(v.1.0.1), ellipsis(v.0.3.2), dplyr(v.1.1.0), backports(v.1.4.1), annotate(v.1.76.0), aod(v.1.3.2), deldir(v.1.0-6), biomaRt(v.2.54.0), vctrs(v.0.5.2), here(v.1.0.1), remotes(v.2.4.2), cachem(v.1.0.6), withr(v.2.5.0), ggforce(v.0.4.1), HDO.db(v.0.99.1), GenomicAlignments(v.1.34.0), treeio(v.1.22.0), prettyunits(v.1.1.1), DOSE(v.3.24.2), ape(v.5.6-2), lazyeval(v.0.2.2), crayon(v.1.5.2), genefilter(v.1.80.3), edgeR(v.3.40.2), pkgconfig(v.2.0.3), tweenr(v.2.0.2), nlme(v.3.1-162), pkgload(v.1.3.2), devtools(v.2.4.5), rlang(v.1.0.6), lifecycle(v.1.0.3), miniUI(v.0.1.1.1), downloader(v.0.4), filelock(v.1.0.2), BiocFileCache(v.2.6.0), rprojroot(v.2.0.3), polyclip(v.1.10-4), graph(v.1.76.0), Matrix(v.1.5-3), aplot(v.0.1.9), boot(v.1.3-28.1), processx(v.3.8.0), png(v.0.1-8), viridisLite(v.0.4.1), rjson(v.0.2.21), bitops(v.1.0-7), gson(v.0.0.9), KernSmooth(v.2.23-20), pander(v.0.6.5), Biostrings(v.2.66.0), blob(v.1.2.3), stringr(v.1.5.0), qvalue(v.2.30.0), remaCor(v.0.0.11), gridGraphics(v.0.5-1), scales(v.1.2.1), memoise(v.2.0.1), GSEABase(v.1.60.0), magrittr(v.2.0.3), plyr(v.1.8.8), gplots(v.3.1.3), zlibbioc(v.1.44.0), compiler(v.4.2.0), scatterpie(v.0.1.8), BiocIO(v.1.8.0), RColorBrewer(v.1.1-3), lme4(v.1.1-31), Rsamtools(v.2.14.0), cli(v.3.6.0), XVector(v.0.38.0), urlchecker(v.1.0.1), patchwork(v.1.1.2), ps(v.1.7.2), MASS(v.7.3-58.2), mgcv(v.1.8-41), tidyselect(v.1.2.0), stringi(v.1.7.12), yaml(v.2.3.7), GOSemSim(v.2.24.0), locfit(v.1.5-9.7), ggrepel(v.0.9.3), grid(v.4.2.0), sass(v.0.4.5), fastmatch(v.1.1-3), tools(v.4.2.0), parallel(v.4.2.0), rstudioapi(v.0.14), foreach(v.1.5.2), gridExtra(v.2.3), farver(v.2.1.1), ggraph(v.2.1.0), digest(v.0.6.31), shiny(v.1.7.4), Rcpp(v.1.0.10), broom(v.1.0.3), later(v.1.3.0), httr(v.1.4.4), AnnotationDbi(v.1.60.0), Rdpack(v.2.4), colorspace(v.2.1-0), brio(v.1.1.3), XML(v.3.99-0.13), fs(v.1.6.1), splines(v.4.2.0), yulab.utils(v.0.0.6), spatstat.utils(v.3.0-1), PROPER(v.1.30.0), tidytree(v.0.4.2), graphlayouts(v.0.8.4), ggplotify(v.0.1.0), plotly(v.4.10.1), sessioninfo(v.1.2.2), xtable(v.1.8-4), jsonlite(v.1.8.4), nloptr(v.2.0.3), ggtree(v.3.6.2), tidygraph(v.1.2.3), ggfun(v.0.0.9), R6(v.2.5.1), RUnit(v.0.4.32), profvis(v.0.3.7), pillar(v.1.8.1), htmltools(v.0.5.4), mime(v.0.12), glue(v.1.6.2), fastmap(v.1.1.0), minqa(v.1.2.5), clusterProfiler(v.4.6.0), BiocParallel(v.1.32.5), codetools(v.0.2-19), fgsea(v.1.24.0), pkgbuild(v.1.4.0), mvtnorm(v.1.1-3), utf8(v.1.2.3), lattice(v.0.20-45), bslib(v.0.4.2), tibble(v.3.1.8), sva(v.3.46.0), pbkrtest(v.0.5.2), curl(v.5.0.0), gtools(v.3.9.4), GO.db(v.3.16.0), survival(v.3.5-3), limma(v.3.54.1), rmarkdown(v.2.20), desc(v.1.4.2), munsell(v.0.5.0), GenomeInfoDbData(v.1.2.9), iterators(v.1.0.14), variancePartition(v.1.28.4), reshape2(v.1.4.4), gtable(v.0.3.1) and rbibutils(v.2.2.13)

message(paste0("This is hpgltools commit: ", get_git_commit()))

## If you wish to reproduce this exact build of hpgltools, invoke the following:

## > git clone http://github.com/abelew/hpgltools.git

## > git reset 465cd7740fb559303acb139de1d478adfc3f6612

## This is hpgltools commit: Fri Mar 24 14:37:21 2023 -0400: 465cd7740fb559303acb139de1d478adfc3f6612

this_save <- paste0(gsub(pattern="\\.Rmd", replace="", x=rmd_file), "-v", ver, ".rda.xz")
message(paste0("Saving to ", this_save))

## Saving to index-v20230330.rda.xz

tmp <- sm(saveme(filename=this_save))