Metadata-Version: 2.1

Name: LeafletSC

Summary: Alternative splicing quantification in single cells with Leaflet

Home-page: https://github.com/daklab/Leaflet

Author: Karin Isaev, Columbia University and NYGC

Requires-Python: >=3.9.15

Description-Content-Type: text/markdown

License-File: LICENSE

Requires-Dist: pyranges==0.0.129

# LeafletSC

LeafletSC is a binomial mixture model designed for the analysis of alternative splicing events in single-cell RNA sequencing data. The model facilitates understanding and quantifying splicing variability at the single-cell level. Below is the graphical model representation:

<p align="center">

</p>

## Compatibility with sequencing platforms

## Getting Started

```bash

```

```Snakemake

{params.regtools_path} junctions extract -a 6 -m 50 -M 500000 {input.bam_use} -o {output.juncs} -s XS -b {output.barcodes}

# Combine junctions and cell barcodes

paste --delimiters='\t' {output.juncs} {output.barcodes} > {output.juncswbarcodes}

```

## Capabilities

With LeafletSC, you can:

- Conduct differential splicing analysis between specific cell groups if cell identities are known.

- Generate synthetic alternative splicing datasets for robust analysis testing.

## If you use Leaflet, please cite our [paper](https://www.biorxiv.org/content/10.1101/2023.10.17.562774v3)

```

2. Add 10X/split-seq mode in addition to smart-seq2

3. Extend framework to seurat/scanpy anndata objects

4. Add notes on generative model and inference method

LeafletSC.egg-info/dependency_links.txt

LeafletSC.egg-info/requires.txt

LeafletSC.egg-info/top_level.txt

LeafletSC/clustering/__init__.py

LeafletSC/clustering/load_cluster_data.py

LeafletSC/clustering/obtain_intron_clusters.py

pyranges==0.0.129

default="no",

help='yes if want to remove lowly used junctions in clusters, default is no')

def process_gtf(gtf_file): #make this into a seperate script that processes the gtf file into gr object that can be used in the main scriptas input

print("The gtf file you provided is " + gtf_file)

print("This step may take a while depending on the size of your gtf file")

# calculate how long it takes to read gtf_file and report it

gtf_exons_gr = gtf_exons_gr.drop_duplicate_positions(strand=True) # Why are so many gone after this?

# Print the number of unique exons, transcript ids, and gene ids

print("The number of unique exons is " + str(len(gtf_exons_gr.exon_id.unique())))

print("The number of unique transcript ids is " + str(len(gtf_exons_gr.transcript_id.unique())))

print("The number of unique gene ids is " + str(len(gtf_exons_gr.gene_id.unique())))

return(gtf_exons_gr)

def filter_junctions_by_shared_splice_sites(df):

#1. Check format of junc_files and convert to list if necessary

# Can either be a list of folders with junction files or a single folder with junction files

# 3. Process each path

for junc_path in junc_files:

junc_path = Path(junc_path)

print(f"Reading in junction files from {junc_path}")

print(f"No junction files found in {junc_path} with suffix {junc_suffix}")

continue

print(f"The number of regtools junction files to be processed is {len(junc_files_in_path)}")

files_not_read = []

try:

juncs = pd.read_csv(junc_file, sep="\t", header=None)

juncs['file_name'] = junc_file # Add the file name as a new column

juncs['cell_type'] = junc_file

all_juncs_list.append(juncs) # Append the DataFrame to the list

except Exception as e:

print(f"Could not read in {junc_file}: {e}")

files_not_read.append(junc_file)

# 5. Concatenate all DataFrames into a single DataFrame

all_juncs = pd.concat(all_juncs_list, ignore_index=True) if all_juncs_list else pd.DataFrame()

all_juncs["intron_length"] = all_juncs["chromEnd"] - all_juncs["chromStart"]

mask = (all_juncs["intron_length"] >= min_intron) & (all_juncs["intron_length"] <= max_intron)

all_juncs = all_juncs[mask]

# Filter for 'chrom' column to handle "chr" prefix

all_juncs = all_juncs.copy()

all_juncs['junction_id'] = all_juncs['chrom'] + '_' + all_juncs['chromStart'].astype(str) + '_' + all_juncs['chromEnd'].astype(str)

# Get total score for each junction and merge with all_juncs with new column "total_counts"

all_juncs = all_juncs.groupby('junction_id').agg({'score': 'sum'}).reset_index().merge(all_juncs, on='junction_id', how='left')

# rename score_x and score_y to total_junc_counts and score

# 9. if singleton is False, remove clusters with only one junction

if singleton == False:

clusters = clusters[clusters.Count > 1]

print("The number of clusters after removing singletons is " + str(len(clusters.Cluster.unique())))

# check if any clusters are singletons now and remove if have singleton == False

if singleton == False:

filtered_clusters_df = filtered_clusters_df[filtered_clusters_df.Count > 1]

print("The number of clusters after removing singletons is " + str(len(filtered_clusters_df.Cluster.unique())))

singleton=True

else:

singleton=False

# print out all user defined arguments that were chosen

print("The following arguments were chosen:" , flush=True)

print("junc_suffix: " + junc_suffix, flush=True)

print("min_num_cells_wjunc: " + str(min_num_cells_wjunc), flush=True)

print("singleton: " + str(singleton), flush=True)

print("filter_low_juncratios_inclust: " + (filter_low_juncratios_inclust), flush=True)

from tqdm import tqdm

import concurrent.futures

import time

parser = argparse.ArgumentParser(description='Read in file that lists junctions for all samples, one file per line and no header')

parser.add_argument('--intron_clusters', dest='intron_clusters',

help='path to the file that has the intron cluster events and junction information from running intron_clustering.py')

help='how you want to name the output file, this will be the input for all Leaflet models')

parser.add_argument('--has_genes', dest='has_genes',

help='yes if intron clustering was done with a gtf file, No if intron clustering was done in an annotation free manner')

help='how many lines to read in at a time, default is 5000')

parser.add_argument('--metadata', dest='metadata',

default=None,

help='path to the metadata file, if provided, the output file will have cell type information')

print("The number of total cells evaluated is " + str(len(all_cells)))

cells_types = clusts[["cell_type", "cell_id"]].drop_duplicates()

print("The number of cells per cell type is:")

print(cells_types.groupby(["cell_type"])["cell_type"].count())

summarized_data = summarized_data.drop_duplicates(subset=['cell_id', 'junction_id'], keep='last') #double check if this is still necessary

summarized_data = clust_cell_counts.merge(summarized_data)

print(np.unique(summarized_data['cell_id'].values))

summarized_data["junc_ratio"] = summarized_data["junc_count"] / summarized_data["Cluster_Counts"]

# if "/" detected in name (cell_type) replace it with "_"

if "/" in name:

name = name.replace("/", "_")

group.to_hdf(output_file + "_" + name + ".h5", key='df', mode='w', complevel=9, complib='zlib')

if metadata is None:

# save summarized_data as hdf file

summarized_data.to_hdf(output_file + ".h5", key='df', mode='w', complevel=9, complib='zlib')

print("Done generating input file for Leaflet model. This process took " + str(round(time.time() - start_time)) + " seconds")

if __name__ == '__main__':

intron_clusters=args.intron_clusters

output_file=args.output_file

has_genes=args.has_genes