# Quick Start

# Quick Start for gimap

For more background on gimap and the calculations done here, read here

    #>
    #> Attaching package: 'dplyr'
    #> The following objects are masked from 'package:stats':
    #>
    #>     filter, lag
    #> The following objects are masked from 'package:base':
    #>
    #>     intersect, setdiff, setequal, union

## Set Up

First we can create a folder we will save files to.

    output_dir <- "output_timepoints"
    dir.create(output_dir, showWarnings = FALSE)

    example_data <- get_example_data("count")

## Setting up data

We’re going to set up three datasets that we will provide to the `set_up()` function to create a `gimap` dataset object.

- `counts` - the counts generated from pgPEN
- `pg_ids` - the IDs that correspond to the rows of the counts and specify the construct
- `sample_metadata` - metadata that describes the columns of the counts including their timepoints

    counts <- example_data %>%
      select(c("Day00_RepA", "Day05_RepA", "Day22_RepA", "Day22_RepB", "Day22_RepC")) %>%
      as.matrix()

`pg_id` are just the unique IDs listed in the same order/sorted the same way as the count data.

    pg_ids <- example_data %>%
      dplyr::select("id")

Sample metadata is the information that describes the samples and is sorted the same order as the columns in the count data.

    sample_metadata <- data.frame(
      col_names = c("Day00_RepA", "Day05_RepA", "Day22_RepA", "Day22_RepB", "Day22_RepC"),
      day = as.numeric(c("0", "5", "22", "22", "22")),
      rep = as.factor(c("RepA", "RepA", "RepA", "RepB", "RepC"))
    )

We’ll need to provide `example_counts`, `pg_ids` and `sample_metadata` to `setup_data()`.

    gimap_dataset <- setup_data(
      counts = counts,
      pg_ids = pg_ids,
      sample_metadata = sample_metadata
    )

It’s ideal to run quality checks first. The `run_qc()` function will create a report we can look at to assess this.

    run_qc(gimap_dataset,
           output_file = "example_qc_report.Rmd",
           overwrite = TRUE,
           quiet = TRUE)

You can take a look at an example QC report here.

    gimap_dataset <- gimap_dataset %>%
      gimap_filter() %>%
      gimap_annotate(cell_line = "HELA") %>%
      gimap_normalize(
        timepoints = "day"
      ) %>%
      calc_gi()

## Example output

Genetic interaction is calculated by:

- `rep` - indicates which sample from the original the data is from. Note the pretreatment is used for calculation and its statistics are not reported.
- `pgRNA_target` - what gene(s) were targeted by this the original pgRNAs for these data
- `mean_expected_cs` - the average expected genetic interaction score
- `mean_gi_score` - the average observer genetic interaction score
- `target_type` - describes whether the CRISPR design is targeting two genes (“gene_gene”), or a gene and a non targeting control (“gene_ctrl”) or a targeting control and a gene (“ctrl_gene”).
- `p_val` - p values from the testing whether a double knockout construct is significantly different in its genetic interaction score from single targets.
- `fdr` - False discovery rate corrected p values

    gimap_dataset$gi_scores %>%
      dplyr::arrange(fdr) %>%
      head() %>%
      knitr::kable(format = "html")

## Plot the results

You can remove any samples from these plots by altering the `reps_to_drop` argument.

    plot_exp_v_obs_scatter(gimap_dataset, reps_to_drop = "Day05_RepA_early")

    plot_rank_scatter(gimap_dataset, reps_to_drop = "Day05_RepA_early")

    plot_volcano(gimap_dataset, reps_to_drop = "Day05_RepA_early", facet_rep = FALSE)

Here’s how you can save plots like the above.

    ggplot2::ggsave("volcano_plot.png")

### Plot specific target pair

We can pick out a specific pair to plot.

    # "NDEL1_NDE1" is top result so let's plot that
    plot_targets_bar(gimap_dataset, target1 = "NDEL1", target2 = "NDE1")

## Saving data to a file

We can save all these data as an RDS or the genetic interaction scores themselves to a tsv file.

    saveRDS(gimap_dataset, "gimap_dataset_final.RDS")

    readr::write_tsv(gimap_dataset$gi_scores, "gi_scores.tsv")

## Session Info

This is just for provenance purposes.

    sessionInfo()