<jupyter_start><jupyter_text>CD4 Cropseq Setup- Hyperparameter check- Select genes to check correlation against- Save the preprocessed AnnData object<jupyter_code>import pandas as pd
import matplotlib.pyplot as plt
import scanpy.api as sc
import scipy as sp
import itertools
import numpy as np
import scipy.stats as stats
from scipy.integrate import dblquad
import seaborn as sns
from statsmodels.stats.multitest import fdrcorrection
import imp
pd.options.display.max_rows = 999
pd.set_option('display.max_colwidth', -1)
import pickle as pkl
import time
import sys
sys.path.append('/wynton/group/ye/mincheol/Github/scrna-parameter-estimation/simplesc')
import scme, utils
data_path = '/wynton/group/ye/mincheol/parameter_estimation/'
with open(data_path + 'interferon_data/immune_genes.pkl', 'rb') as f:
    immune_genes = pkl.load(f)
with open(data_path + 'cd4_cropseq_data/genes_to_test.pkl', 'rb') as f:
    genes_to_test = pkl.load(f)
genes_to_test<jupyter_output><empty_output><jupyter_text>Read the AnnData object<jupyter_code>adata = sc.read(data_path + 'cd4_cropseq_data/guide_singlets.h5ad')<jupyter_output><empty_output><jupyter_text>For each group, compare the KO gene with the downstream genes<jupyter_code>test_guides = ['WT', 'FUBP1.77964146', 'ZNF706.101201657']
adata = adata[adata.obs.group.isin(test_guides)]
adata.shape
estimator = scme.SingleCellEstimator(
    adata=adata, 
    group_label='group',
    n_umis_column='n_counts',
    num_permute=10000,
    beta=0.1)
estimator.estimate_beta_sq(tolerance=3)
%matplotlib inline
estimator.plot_cv_mean_curve()
for group in adata.groups.drop_duplicates()[:2]:
    group
['WT', 'FUBP1.77964146', 'ZNF706.101201657', 'RFX2.6026185', 'VDR.47855738']
adata.obs.group.value_counts().index[:5]
adata.obs.head(5)<jupyter_output><empty_output><jupyter_text>Select relevant cells and save the AnnData object<jupyter_code>adata = sc.read('/wynton/group/ye/rgate/072319_tfko_140_seq/tfko.sng.guides.h5ad')
df = pd.read_csv(
    '/wynton/group/ye/rgate/072319_tfko_140_seq/tfko.sng.guides.combat.n50.res1.meta.txt', 
    sep='\t',
    index_col=0)
df_singlet = df.query('guide_multiplet_cov == "1"').copy()
df_singlet['group'] = df_singlet.apply(lambda row: '.'.join(row['guide1_cov'].split('.')[:-1]) if row['WT'] == 'F' else 'WT', axis=1)
adata_singlets = adata[df_singlet.index].copy()
adata_singlets.obs['group'] = df_singlet['group']
ko_genes = list(set([x.split('.')[0] for x in adata_singlets.obs.guide1_cov.drop_duplicates().tolist()]))
genes = list(set(ko_genes) | set(immune_genes))
adata_singlets = adata_singlets[:, genes].copy()
adata_singlets.shape
sc.write(data_path + 'cd4_cropseq_data/guide_singlets.h5ad', adata_singlets)<jupyter_output>... storing 'group' as categorical<jupyter_text>Get a list of genes to test correlations against<jupyter_code>prefixes = ['STAT', 'IRF', 'IFI', 'ETS', 'NCOA', 'JUN', 'BATF', 'SATB1', 'FOXP']
genes_to_test = []
for gene in ko_genes:
    for prefix in prefixes:
        if prefix in gene:
            genes_to_test.append(gene)
with open(data_path + 'cd4_cropseq_data/genes_to_test.pkl', 'wb') as f:
    pkl.dump(genes_to_test, f)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.13703820650457846, 'fraction_numerical': 0.03473318598042311, 'mean_word_length': 6.167420814479638, 'pattern_counts': {'":': 0, '<': 20, '<?xml version=': 0, '>': 20, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '5', 'dir': 'jupyter-structured-clean-dedup', 'id': '179143', 'license': 'MIT', 'n_tokens_mistral': 1317, 'n_tokens_neox': 1224, 'n_words': 223, 'path': 'analysis/cd4_cropseq/cd4_cropseq_setup.ipynb', 'repo_name': 'yelabucsf/scrna-parameter-estimation', 'size': '3167'}