FlowCode_Unmixing_Workflow.Rmd

---
title: "FlowCode Spectral Flow Unmixing Workflow"
output: html_notebook
---

This is an [R Markdown](http://rmarkdown.rstudio.com) Notebook. When you execute code within the notebook, the results appear beneath the code. 

Try executing this chunk by clicking the *Run* button within the chunk or by placing your cursor inside it and pressing *Ctrl+Shift+Enter*. 

```{r}
2 + 2
```

When you save the notebook, an HTML file containing the code and output will be saved alongside it (click the *Preview* button or press *Ctrl+Shift+K* to preview the HTML file).

The preview shows you a rendered HTML copy of the contents of the editor. Consequently, unlike *Knit*, *Preview* does not run any R code chunks. Instead, the output of the chunk when it was last run in the editor is displayed.

# Installation

You will both need `AutoSpectral` and `FlowCodeUnmix` to run this, plus their compiled accelerators, `AutoSpectralRcpp` and `FlowCodeUnmixRcpp`. The accelerators are no longer just a speed convenience: `get.flowcode.spectra()` now *requires* `AutoSpectralRcpp` (there's no pure-R fallback for it), and the staged pipeline used by `unmix.flowcode.fcs()`/`unmix.flowcode.folder()` requires `AutoSpectralRcpp >= 1.1.0` specifically.

The easiest way to get everything is [`pak`](https://pak.r-lib.org/). `pak::pkg_install()` resolves CRAN, Bioconductor (via the `bioc::` prefix) and GitHub (`"owner/repo"`) packages together in one call, with proper dependency resolution and parallel downloads -- there's no separate `BiocManager::install()` or `remotes::install_github()` step needed.

```{r}
# one-time, if you don't already have it
install.packages("pak")

pak::pkg_install(c(
  "bioc::flowWorkspace",
  "bioc::FlowSOM",
  "shiny",
  "shinyFiles",
  "ggplot2",
  "scattermore",
  "dplyr",
  "DrCytometer/AutoSpectral",
  "DrCytometer/FlowCodeUnmix",
  "DrCytometer/AutoSpectralRcpp",
  "DrCytometer/FlowCodeUnmixRcpp"
))
```

Once everything is installed, the Rcpp accelerators are picked up automatically wherever they help, no further action needed.


# Basic AutoSpectral workflow

We start as we would with any spectral flow cytometry data set in AutoSpectral, using the automated spectral extraction pipeline. We will load the cytometer-specific parameters, create and check the control file, extract the fluorophore spectra directly from the raw single-stained controls, identify autofluorescence signatures and measure the spectral variation in the fluorophore output. For the full walkthrough of this pipeline, including all of the diagnostic plots it produces, see [AutoSpectral Workflow](https://drcytometer.github.io/AutoSpectral/articles/01_Automated_Spectral_Workflow.html).


## Load the parameters

Here we are using data from the Cytek Aurora (5-laser system). For data from a different spectral cytometer, call the appropriate cytometer as "cytometer". For more details on the AutoSpectral workflow, see [AutoSpectral help](https://drcytometer.github.io/AutoSpectral/).

```{r, startup}
library(AutoSpectral)
autospectral.parameters <- AutoSpectral::get.autospectral.param(
  cytometer = "aurora",
  figures = TRUE
)
```


## Create and check the control file

The control file is a CSV file (spreadsheet) that describes your control files for AutoSpectral. A draft of this will be created automatically by AutoSpectral if you tell it where your single-stained control files are. You will need to edit this and check it to be sure it is correct. Details on how to do this are in this [article](https://drcytometer.github.io/AutoSpectral/articles/02_Control_File_example.html). There is also a [Shiny helper](https://github.com/DrCytometer/AutoSpectralHelper) application for this, which opens as an interactive webpage.

```{r, define the controls}
# where are your single-stained control files?
# replace the "my_control_folder" with the path to your files.
# Tip: if you start RStudio by double-clicking on this or another R document,
# that tells RStudio where you want to work. If you place a folder (say, called SSC)
# containing your control fcs files in that folder, you will be able to access 
# it with just "./SSC
control.folder <- "C://User/Documents/my_control_folder"

# or
control.folder <- "./SSC"

AutoSpectral::create.control.file(
  control.dir = control.folder,
  asp = autospectral.parameters
)
```

`create.control.file()` now defaults to `legacy = FALSE`, so the draft control file is set up for the automated extraction pipeline below: you'll still need to fill in `marker`, `fluorophore`, `control.type`, and `universal.negative`, but you can leave `gate.name`/`gate.define` blank -- those columns are only read by the legacy gating pipeline. The peak-channel `channel` column is optional too; it's determined automatically from the data.

Once you have edited your control file, check it:

```{r, check the controls}
control.file <- "fcs_control_file.csv"
AutoSpectral::check.control.file(
  control.dir = control.folder,
  control.def.file = control.file,
  asp = autospectral.parameters
)
```


## Extract the fluorophore spectral profiles

This replaces the old three-step `define.flow.control()` -> `clean.controls()` -> `get.fluorophore.spectra()` sequence with a single call. No manual gating is required: saturation removal, singlet gating, AF-orthogonalised cosine filtering of the brightest candidate events, and kNN scatter-matched per-event background subtraction all happen internally, with an automatic per-fluorophore fallback to the legacy pipeline if the automated result doesn't pass its own cosine-similarity QC check against the reference spectral library.

```{r, get the spectral profiles}
fluorophore.spectra <- AutoSpectral::get.spectra.automated(
  control.dir = control.folder,
  control.def.file = control.file,
  asp = autospectral.parameters
)
```

The default settings work well across a wide range of panels, but can be tuned if needed (`n.candidates`, `n.spectral`, `k.neighbors`, `cosine.threshold`, etc. -- see [AutoSpectral Workflow](https://drcytometer.github.io/AutoSpectral/articles/01_Automated_Spectral_Workflow.html) for the full parameter reference). Check the diagnostic plots in `figure_spectra`, `figure_similarity_heatmap`, `figure_scatter_match` and the per-fluorophore QC report in `figure_spectral_ribbon`; anything flagged `LEGACY_REFINEMENT` in the console QC summary is worth a look.

We also need a minimal `flow.control` object -- this just carries the antigen/fluorophore name lookup used later on (e.g. by `unmix.backbone()`, `launch.threshold.app()`, and when writing FCS keywords). It replaces the old, fully-gated object that used to come out of `define.flow.control()`/`clean.controls()`.

```{r, minimal flow control}
flow.control <- AutoSpectral::reload.flow.control(
  control.dir = control.folder,
  control.def.file = control.file,
  asp = autospectral.parameters
)
```


## Measure autofluorescence spectra

Autofluorescence is always present in our stained cells, including our single-stained cell controls. This is a source of variability that contributes to noise (spread) in our unmixed data. We need to identify the variability in autofluorescence to counteract this. With FlowCode data, we will often have different tissue sources, each of which will have a distinct distribution of autofluorescence profiles. We will eventually need to identify all of those. To start, though, we will identify only the variability in the type of sample we are using for the single-stained controls. For instance, with mouse studies, these will likely be splenocytes, so we will use unstained spleen as the source of our autofluorescence.

In this example, I have an unstained spleen sample as the "Unstained" (so named by SpectroFlo) in my reference control group. This FCS file is in the single-stained control folder, along with all of the other controls.

```{r, measure AF}
spleen.autofluorescence <- AutoSpectral::get.af.spectra(
  unstained.sample = file.path( control.folder, "Unstained.fcs" ),
  asp = autospectral.parameters,
  spectra = fluorophore.spectra
)
```

Note that `refine` currently defaults to `TRUE`: after the initial SOM-based clustering, a second pass of targeted modulation runs on cells that remain furthest from zero, which helps with messier tissue samples at the cost of some extra time. If you're working with clean lymphoid samples (spleen, PBMCs) and want the quicker first-pass-only result, set `refine = FALSE`.


## Determine spectral variation in the fluorophores

Now that we know the contribution of autofluorescence to the variation, we can extract the variation from the fluorophore emissions themselves.

```{r, measure variation}
fluorophore.variation <- AutoSpectral::get.spectral.variants(
  control.dir = control.folder,
  control.def.file = control.file,
  asp = autospectral.parameters,
  spectra = fluorophore.spectra,
  parallel = FALSE
)
```

As of AutoSpectral 1.6.0, `get.spectral.variants()` no longer takes an `af.spectra` argument -- autofluorescence for this step is now derived in situ, per fluorophore, from the universal-negative files already listed in your control table. Passing `af.spectra` here (even positionally) will trigger a deprecation warning and be silently ignored, so don't rely on it.

That is the end of the normal AutoSpectral workflow for us. We now need to do some FlowCode-specific steps before proceeding to the debarcoding and unmixing.


# FlowCode Unmixing Workflow

For the FlowCode workflow, we need some additional pieces of information:
* Which are the FlowCode epitope tag channels?
* What are the valid FlowCode combinations (in the combination file)?
* What is the name (guide/TCR) of each combination (in the combination file)?
* Where are the cells (automatically gated around the brightest FlowCode-expressing events)?
* What variation in fluorophore output is created by the combinations due to FRET or a FRET-like artifact? For this, we need your backbone control, stained only with the FlowCode epitope antibodies. All combinations should be well represented for best results.


## Where is everything?

Provide the location (file path) and name of your backbone FCS file. Provide the name (and file path) of your combination CSV file.

```{r, FlowCode combos}
flowcode.backbone <- "./Raw/FC control (Spleen)/F2 FC_FRET_02_Plate_002.fcs" 
combo.file <- "OTB01 Treg Tx factor combinations.csv"
```


Now we can unmix the backbone control, using a slight adaptation of the AutoSpectral unmixing approach. This is a spleen-based sample, so we use `spleen.autofluorescence`.

```{r, load the backbone data}
FlowCodeUnmix::unmix.backbone(
   flowcode.backbone.fcs = flowcode.backbone,
   spectra = fluorophore.spectra,
   af.spectra = spleen.autofluorescence,
   spectra.variants = fluorophore.variation,
   flow.control = flow.control,
   asp = autospectral.parameters,
   flowcode.combo.file = combo.file,
   output.dir = "./flowcode_spectra",
   filename = "FlowCode_Backbone.rds",
   parallel = TRUE
)
```

A couple of things changed here beyond the `flow.control` rename above:

* `spectra.variants` is a new, optional-but-recommended argument. When supplied, each combo fluorophore's per-cell best-fitting spectral variant is substituted before the FRET residual is computed (identically for the combo's own cells and for the untransduced background pool), so ordinary single-colour spectral drift isn't mistaken for FRET. It's also a prerequisite, together with `af.spectra`, for `plot.corrections = TRUE` later in `get.flowcode.spectra()`.
* Before the gate is even defined, the backbone data now goes through automatic saturation removal and a two-pass Area/Height scatter-ratio singlet gate (`singlet.quantiles = c(0.85, 0.975)` by default, skipped with a message if Height channels aren't present). A doublet carrying two different barcodes would otherwise still debarcode to a single Id, and its unmixing residual -- driven by the second, uncounted barcode's own spillover -- would look exactly like FRET and contaminate the measured FRET spectra.
* `parallel`/`threads` now control the AF-extraction unmix over the gated events, same convention as elsewhere (`threads = 0` means "all available cores").

This saves the output as an RDS file (a large one, plus a downsampled `Small_...` companion capped at 30,000 events) in the folder specified by `output.dir`. If you change `output.dir`, you'll need to change it later in `get.flowcode.spectra()`.
*Note: change to returning object

## Setting thresholds

This part requires you to manually select the point on the graph where the data start becoming "positive" for each of the FlowCode epitope tags. That is, where is the threshold for having the FlowCode or not? This is manual because it's pretty hard to do it well in an automated manner considering that at this stage we still have FRET-based unmixing issues.

I have created a Shiny app to allow you to do this interactively, as in Orian Bricard's FlowCode debarcoding app.

```{r, run app}
FlowCodeUnmix::launch.threshold.app(
  backbone.rds = "./flowcode_spectra/FlowCode_Backbone.rds"
)
```

Pointing `backbone.rds` at the RDS file saved by `unmix.backbone()` lets the app load the fluorophore-tag mapping and spectra straight out of it, so you no longer need to separately pass `flowcode.combo.file`/`flow.control` here (the downsampled `Small_...` file works just as well and loads faster). If you'd rather browse for a file interactively, or point the app at a raw/not-yet-unmixed FCS file, you can still leave `backbone.rds` as `NULL` and use the app's own file picker.

When you're done, you should have a CSV file (spreadsheet) with `Fluor` and `Threshold_Raw` columns, giving the threshold for each FlowCode channel in raw, untransformed unmixed scale.


## Measuring FRET (unmixing) errors due to FlowCodes

In this step, we will take the FlowCode epitope-stained backbone data (from the unmixed data in the RDS file), debarcode it to identify valid combinations, and assess spectral unmixing errors per combination. This gives us a set of potential corrections to apply to each cell for a given barcode combination.

```{r, find FRET errors}
flowcode.spectra <- FlowCodeUnmix::get.flowcode.spectra(
  backbone.rds = "./flowcode_spectra/FlowCode_Backbone.rds",
  thresholds.file = "./flowcode_spectra/thresholds.csv",
  asp = autospectral.parameters,
  spectra.variants = fluorophore.variation,
  af.spectra = spleen.autofluorescence,
  output.dir = "./flowcode_spectra",
  filename = "FlowCode_Spectra.rds",
  plot.corrections = TRUE
)
```

This function has changed substantially. Note in particular:

* **`AutoSpectralRcpp` is now a hard requirement** -- Stage 1 (a debarcoding-quality unmix, using the same tuning the production pipeline uses so characterisation and correction agree on debarcoding/background/AF from the start) runs through `AutoSpectralRcpp::unmix.autospectral.rcpp()`. There is currently no pure-R fallback for this function.
* `spectra.variants` and `af.spectra` are now required arguments to get `plot.corrections = TRUE` to work, since the before/after plots re-run the full production pipeline (`unmix.flowcode.cpp.staged()`) on the backbone.
* Background is now handled by a single **pooled, cross-combo** median correction per FlowCode fluorophore (computed from off-target observations across every valid combo, then subtracted from every cell, transduced or not), rather than anything per-combo -- this is what separates genuine non-specific background/staining from real FRET.
* `background.pool.n` (default `5000`) controls how many untransduced cells are used as the scatter-matching/reference pool, and `n.cells.per.combo` (default `500`) caps how many of a combo's own brightest events are used when measuring its FRET spectra.
* `som.dim` (default `10`) is the base SOM grid side length used to cluster each combo's residuals into candidate FRET variants; it's shrunk automatically for combos with fewer than 500 retained cells.

There are also two optional QC filters, both off by default, that you can turn on if you see noisy or implausible FRET variants coming out of a combo:

```{r, find FRET errors with QC}
flowcode.spectra <- FlowCodeUnmix::get.flowcode.spectra(
  backbone.rds = "./flowcode_spectra/FlowCode_Backbone.rds",
  thresholds.file = "./flowcode_spectra/thresholds.csv",
  asp = autospectral.parameters,
  spectra.variants = fluorophore.variation,
  af.spectra = spleen.autofluorescence,
  output.dir = "./flowcode_spectra",
  plot.corrections = TRUE,
  fret.median.qc = TRUE,               # loose sanity filter vs. the combo's own median FRET direction
  fret.median.sim.threshold = 0,       # 0 = only excludes events pointing the opposite way
  fret.background.qc = TRUE,           # screens candidate variants against background-pool contamination
  fret.background.sim.threshold = 0.8  # absolute cosine similarity cutoff for that screen
)
```

`fret.median.qc` drops events, before SOM clustering, that point away from the combo's own median FRET direction -- it's a loose filter, not a same-shape one, since a combo's true FRET isn't guaranteed to have a single direction. `fret.background.qc` instead screens the SOM-derived candidate variants (never the median row) after clustering, against the FRET-shaped signal found in the untransduced background pool run through the same cleaning; variants that closely resemble that background contamination are dropped. They address different failure modes and can be used together.

This call saves the output as an RDS file in the folder specified by `output.dir`. If you select `plot.corrections = TRUE`, it will identify the FRET errors on a cell-by-cell basis, correct them, and plot examples of what the FlowCode epitope unmixing looks like on the backbone file with and without corrections. That can be a bit slow, but is definitely worth doing the first couple of times.


## FlowCode Unmixing

We now have all the data required for unmixing, or, at least we do for the any files coming from spleen. We'll get to other tissue sources in a second.

To unmix a single file, call `unmix.flowcode.fcs()`.

```{r, unmixing}
# for example:
fully.stained.spleen <- "./Raw/Spleen/Spleen_Mouse_001.fcs"

# then call:
FlowCodeUnmix::unmix.flowcode.fcs(
  fcs.file = fully.stained.spleen,
  spectra = fluorophore.spectra,
  asp = autospectral.parameters,
  flow.control = flow.control,
  af.spectra = spleen.autofluorescence,
  spectra.variants = fluorophore.variation,
  flowcode.spectra = flowcode.spectra,
  parallel = TRUE
)
```

By default this runs the new **staged** pipeline (`unmix.flowcode.cpp.staged()`): a generic AutoSpectral unmix for debarcoding, the pooled cross-combo background subtraction described above, a small joint FRET fit against the background-corrected raw data, and a final generic AutoSpectral unmix for the reported values. This requires `AutoSpectralRcpp >= 1.1.0`; if it's missing or out of date you'll get an explicit error telling you to update it, and if the C++ call itself fails for some other reason, it falls back with a warning to the (currently outdated) pure-R `unmix.flowcode()`, which doesn't yet support the extra tuning parameters below.

If you want to tune the staged pipeline, the key extra parameters (all with sensible defaults, so you only need to touch the ones you care about) are:

```{r, unmixing with tuning}
FlowCodeUnmix::unmix.flowcode.fcs(
  fcs.file = fully.stained.spleen,
  spectra = fluorophore.spectra,
  asp = autospectral.parameters,
  flow.control = flow.control,
  af.spectra = spleen.autofluorescence,
  spectra.variants = fluorophore.variation,
  flowcode.spectra = flowcode.spectra,
  n.passes = 1,               # rounds of per-cell spectral optimization
  n.af.passes = 1,            # rounds of AF extraction
  refine.af.quantile = 0.5,   # passed through to unmix.autospectral.rcpp()
  alpha = 0.5,
  collinear.threshold = 0.5,
  joint.pair.resolution = TRUE,
  fret.median.only = FALSE,   # TRUE skips the per-cell FRET variant scan, always uses the median FRET row
  return.diagnostics = FALSE, # TRUE returns per-cell diagnostics instead of NULL (see below)
  chunk.size = 2e6,           # events per chunk, for memory management on large files
  parallel = TRUE
)
```

`return.diagnostics = TRUE` is worth knowing about if you're debugging a specific combo or fluorophore: instead of writing the file and returning `NULL`, you get back (invisibly) a list with `median.background`, `flowcode.ids`, `fret.k`, `fret.index`, `resid.ratio` and `leakage.ratio` -- note `median.background` is pooled per chunk rather than across the whole file, so treat it as indicative rather than exact for very large files processed in several chunks.

To unmix a folder containing many files, call `unmix.flowcode.folder()`.

```{r, unmixing a bunch of stuff}
# for example:
spleen.folder <- "./Raw/Spleen/"

# then call:
FlowCodeUnmix::unmix.flowcode.folder(
  fcs.dir = spleen.folder,
  spectra = fluorophore.spectra,
  asp = autospectral.parameters,
  flow.control = flow.control,
  af.spectra = spleen.autofluorescence,
  spectra.variants = fluorophore.variation,
  flowcode.spectra = flowcode.spectra,
  parallel = TRUE
)
```

If we want to look at other tissues sources, or any samples where the autofluorescence may differ from the autofluorescence profiles we have already extracted, we first need to call `AutoSpectral::get.af.spectra()` on unstained samples representing those autofluorescence sources. We can then unmix samples from those tissues, extracting the matching autofluorescence correctly.

```{r, unmixing different tissues}
# Let's say we have brain samples

# this is our unstained brain sample (raw data)
unstained.brain <- "./Raw/Unstained controls/Unstained brain.fcs"

# we extract the AF profiles from it
brain.autofluorescence <- AutoSpectral::get.af.spectra(
  unstained.sample = unstained.brain,
  asp = autospectral.parameters,
  spectra = fluorophore.spectra
)

# now we want the raw data for the fully stained brain sample
fully.stained.brain <- "./Raw/Brain/Brain_Mouse_001.fcs"

# then call to unmix:
FlowCodeUnmix::unmix.flowcode.fcs(
  fcs.file = fully.stained.brain, # change the sample
  spectra = fluorophore.spectra,
  asp = autospectral.parameters,
  flow.control = flow.control,
  af.spectra = brain.autofluorescence, # change the autofluorescence to match
  spectra.variants = fluorophore.variation,
  flowcode.spectra = flowcode.spectra,
  parallel = TRUE
)
```

## Looping over several tissues

Doing the above by hand for every tissue gets old fast once you have more than one or two sources. `unmix.flowcode.folder()` unmixes an entire directory of FCS files against a single `af.spectra`, so the natural setup is **one subfolder per tissue** for your fully-stained files, plus **one unstained FCS file per tissue** to build that tissue's AF spectra from. Something like:

```
Raw/
├── Unstained controls/
│   ├── Unstained spleen.fcs
│   ├── Unstained brain.fcs
│   └── Unstained liver.fcs
├── Spleen/
│   ├── Spleen_Mouse_001.fcs
│   └── Spleen_Mouse_002.fcs
├── Brain/
│   ├── Brain_Mouse_001.fcs
│   └── Brain_Mouse_002.fcs
└── Liver/
    ├── Liver_Mouse_001.fcs
    └── Liver_Mouse_002.fcs
```

The important thing is that a single `unmix.flowcode.folder()` call is all-or-nothing for a directory: every FCS file it finds there gets unmixed with the same `af.spectra`, so don't mix tissues within one folder. It's worth mirroring that structure on the output side too (one output subfolder per tissue), so you're never relying on filename suffixes alone to tell tissues apart later.

With that layout, key the unstained file, the stained folder, and the output folder by the same tissue name, so the loop matches them up by name rather than by position (much safer if you ever add or reorder a tissue):

```{r, tissue lookup tables}
# named vector: tissue name -> unstained FCS file used to build that tissue's AF spectra
tissue.unstained <- c(
  Spleen = "./Raw/Unstained controls/Unstained spleen.fcs",
  Brain  = "./Raw/Unstained controls/Unstained brain.fcs",
  Liver  = "./Raw/Unstained controls/Unstained liver.fcs"
)

# named vector: tissue name -> folder of fully-stained FCS files for that tissue
tissue.folders <- c(
  Spleen = "./Raw/Spleen/",
  Brain  = "./Raw/Brain/",
  Liver  = "./Raw/Liver/"
)

# named vector: tissue name -> output folder for that tissue's unmixed files
tissue.output <- c(
  Spleen = "./Unmixed/Spleen/",
  Brain  = "./Unmixed/Brain/",
  Liver  = "./Unmixed/Liver/"
)
```

Then the loop itself just does, per tissue, what we did by hand above: extract that tissue's AF spectra, then unmix its whole folder against it.

```{r, tissue loop}
for ( tissue in names( tissue.unstained ) ) {

  message( "Processing tissue: ", tissue )

  # 1. extract this tissue's AF profile from its unstained sample
  tissue.af <- AutoSpectral::get.af.spectra(
    unstained.sample = tissue.unstained[[ tissue ]],
    asp = autospectral.parameters,
    spectra = fluorophore.spectra
  )

  # 2. unmix every fully-stained file for this tissue against that AF
  FlowCodeUnmix::unmix.flowcode.folder(
    fcs.dir = tissue.folders[[ tissue ]],
    spectra = fluorophore.spectra,
    asp = autospectral.parameters,
    flow.control = flow.control,
    af.spectra = tissue.af,
    spectra.variants = fluorophore.variation,
    flowcode.spectra = flowcode.spectra,
    output.dir = tissue.output[[ tissue ]],
    parallel = TRUE
  )
}
```

A couple of practical notes on this pattern:

* If you want to keep each tissue's `af.spectra` around afterward (for later re-use, or for debugging a specific tissue), stash it in a list instead of letting the loop variable get overwritten each iteration: `tissue.af.list <- list(); tissue.af.list[[ tissue ]] <- tissue.af`.
* `lapply` over `names( tissue.unstained )` works just as well as the `for` loop if you prefer that style -- there's nothing about this pattern that depends on the loop body running sequentially.
* This uses the same `spectra`, `flow.control`, `spectra.variants` and `flowcode.spectra` for every tissue -- only `af.spectra` and the file locations change. That's deliberate: the fluorophore spectra, spectral variants and the FRET correction library are all properties of your reagents/backbone, not of the tissue, and shouldn't be recomputed per tissue.
* `unmix.flowcode.folder()` doesn't expose the staged-pipeline tuning parameters (`cell.weight`, `noise.floor`, `alpha`, etc. -- see the note above), so if you need those per-tissue too, adapt the loop to call `unmix.flowcode.fcs()` over `list.files( tissue.folders[[ tissue ]], pattern = ".fcs", full.names = TRUE )` instead of `unmix.flowcode.folder()`.
* As the next section covers, positivity thresholds can also drift by tissue -- if you see that, add a `thresholds.file` lookup (keyed by tissue name the same way) alongside `tissue.unstained`/`tissue.folders`/`tissue.output` and pass `thresholds.file = tissue.thresholds[[ tissue ]]` into the unmixing call.

In many cases, the thresholds for positivity may vary slightly by tissue type or sample source. If that is the case, you will need to provide new thresholds as part of the unmixing process. To start, unmix a single FCS file from the source you want to use to set new thresholds. Do the unmixing with AutoSpectral, extracting the AF per cell, since this will get you closer to what the data will look like with the FlowCodeUnmix. Then load this file into the ThresholdApp, select new thresholds for each FlowCode marker, and save the results (you can and should rename the CSV file). Now provide this new CSV file containing the updated thresholds to the unmixing.

I suspect that in most cases there will be little need for re-thresholding between different tissues. There may be a need to re-threshold on a fully stained sample, though, since this will differ from the backbone control we have used to establish the first set of thresholds.

```{r, changing the thresholds}

# this is our unstained brain sample (raw data)
unstained.brain <- "./Raw/Unstained controls/Unstained brain.fcs"

# we extract the AF profiles from it
brain.autofluorescence <- AutoSpectral::get.af.spectra(
  unstained.sample = unstained.brain,
  asp = autospectral.parameters,
  spectra = fluorophore.spectra
)

# now we want the raw data for the fully stained brain sample
fully.stained.brain <- "./Raw/Brain/Brain_Mouse_001.fcs"


FlowCodeUnmix::launch.threshold.app(
  flowcode.combo.file = combo.file,
  flow.control = flow.control,
  spectra = fluorophore.spectra,
  asp = autospectral.parameters
)

new.thresholds <- "./flowcode_spectra/New_thresholds.csv"

# then call to unmix:
FlowCodeUnmix::unmix.flowcode.fcs(
  fcs.file = fully.stained.brain, # change the sample
  spectra = fluorophore.spectra,
  asp = autospectral.parameters,
  flow.control = flow.control,
  af.spectra = brain.autofluorescence, # change the autofluorescence to match
  spectra.variants = fluorophore.variation,
  flowcode.spectra = flowcode.spectra,
  thresholds.file = new.thresholds,
  parallel = TRUE
)

```

Note that `flowcode.combo.file` and `flow.control` must be supplied to `launch.threshold.app()` together, or not at all -- a partial pair now raises an explicit error up front, rather than failing confusingly later on. As above, you can also just point `backbone.rds` at any RDS produced by `unmix.backbone()` and skip supplying these separately.

# FlowCode Debarcoding Workflow

With the new unmixing, the resulting FCS files are already debarcoded. If you inspect the files in R or FlowJo, you should see a channel for each of the tags (names, guides, TCRs) you have linked to each FlowCode combination. The expression in these channels is the expression of that combination of FlowCode epitopes in the cells. Additionally, you should have a "FlowCode" channel, which is the expression level of the epiopes in cells with valid FlowCode combinations. In other words, the "FlowCode" channel measures transduction.

So, in FlowJo, you can gate on your cells and then plot a tag combo versus whatever markers you have in your flow panel. This is useful for checking CRISPR guide efficiency (for example, plotting GATA-3 expression in Gata3-targeted cells) or for phenotyping. If you want, you can run dimensionality reduction and clustering approaches in FlowJo, R or Python on the FCS files as they are or with some pre-gating.

For our analyses, we recommend pre-gating to specific biologically relevant cell populations. I still need to write the part to do this in R. At present, you should be able to treat the unmixed FCS files you get from this workflow as you would previously, running them through Orian Bricard's [FlowCodeDecoder](https://github.com/DrCytometer/FlowcodeDecoder). The difference is that the unmixing will now be correct (hopefully), both for the barcodes and for the phenotyping and gating markers.