<jupyter_start><jupyter_text>16_snp_invocab4/9/21check if vocabl selected regions ar ein any SNPs associated with any disease (ld expansion)<jupyter_code>import os, glob
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import pybedtools
vocab_regions_files = glob.glob('../data/processed/fig4_modelling/vocab_sum_final/vocab_instances_df_filt_*csv')
vocab_regions_files

def get_snp_vocab(snp_bed):
    snp_vocab_all = pd.DataFrame()
    for file in vocab_regions_files:
        tissue = os.path.basename(file).split('_')[-1].split('.')[0]
        print(tissue)
        vocab_regions_df = pd.read_csv(file)
        vocab_regions_df1 = vocab_regions_df[['chr','vocab1_start','vocab1_stop','vocab1','vocab_pair']]
        vocab_regions_df1.columns = ['chr','start','stop','vocab','vocab_pair']
        vocab_regions_df2 = vocab_regions_df[['chr','vocab1_start','vocab1_stop','vocab1','vocab_pair']]
        vocab_regions_df2.columns = ['chr','start','stop','vocab','vocab_pair']
        vocab_regions_df = pd.concat([vocab_regions_df1, vocab_regions_df2])
        vocab_regions_bed = pybedtools.BedTool.from_dataframe(vocab_regions_df)
        intersect_bed = vocab_regions_bed.intersect(snp_bed, wa=True, wb=True)
        print(len(intersect_bed))
        if len(intersect_bed)>0:
            intersect_df = intersect_bed.to_dataframe().drop_duplicates()
            print(intersect_df.shape)
            intersect_df['tissue'] = tissue
            snp_vocab_all = pd.concat([snp_vocab_all, intersect_df])

    snp_vocab_all.columns = ['chr_vocab','start_vocab','stop_vocab','tf_name','vocab_pair','chr_snp','start_snp','stop_snp','Linked_SNP','tissue']
    print(snp_vocab_all.shape)
    return snp_vocab_all<jupyter_output><empty_output><jupyter_text>Try 1: older version of eQTL linked snp filtered by DHS peaks (28 hits)Jan 2019 vs Nov 2018 15 --> 28 hits<jupyter_code># snp_df = pd.read_csv('../data/external/cis_eQTLs_all_re.txt',sep='\t')


snp_df = pd.read_csv('../data/external/eQTL-linked-index-SNP-DHS-disease-IndexP-Jan102019.txt',sep='\t')
snp_bed_df = snp_df[['Chr','Position','Linked_SNP']]
snp_bed_df['stop'] = snp_bed_df.Position+1
snp_bed_df = snp_bed_df[['Chr','Position','stop','Linked_SNP']]
snp_bed_df.columns = ['chr','start','stop','name']



snp_bed = pybedtools.BedTool.from_dataframe(snp_bed_df)
print(snp_df.shape)
display(snp_df[:5])

snp_vocab_all = get_snp_vocab(snp_bed).drop_duplicates()
print(snp_vocab_all.shape)

snp_vocab_all = snp_vocab_all.merge(snp_df[['Linked_SNP','Index_SNP','Disease','DHS_tissue','eGene','eQTL_tissue']]).drop_duplicates()
snp_vocab_all[-40:]#.columns

# # these are the ones we actually tested DO NOT RERUN

# snp_vocab_all[['chr_vocab', 'start_vocab', 'stop_vocab', 'tf_name', 'vocab_pair',
#        'chr_snp', 'start_snp', 'stop_snp', 'Linked_SNP', 'tissue']].drop_duplicates()
snp_vocab_all.to_csv('../data/processed/snp_motifs/snp_vocab_all_wdhs_egene.csv')<jupyter_output><empty_output><jupyter_text>foray into actual stuffthese are the ones we actually tested apparentlychr14	69253364	69253380	MAFG	KLF5::MAFG	chr14	69253364	69253365	rs12434551	GDSD65498	KLF5::MAFG	86	59	chr14	69253364	69253450	69253340	69253474	MAFG	69253364	69253380	KLF5	69253439	69253450	GDSD6	loop_loop	2	KLF5::MAFG+GDSD6+2	0	TRUE	scrA_B_seq	CTTCGGAAGAGGCCTGGTGGAATACGCCTATAGTGATAGCGTCCCCAAAACTTCCTGGGTGTAAAGGAACCAAAGGCAAGTGAGCTGCAAATAGCAATCTCTGCCTCTCCCCAAATAAATCCTCCCAAGGCTCAC	KLF5::MAFG+GDSD6+2+scrA_BMAFG was the one with the snp rs12434551 init (GDSD6)35	chr14	23289048	23289065	PLAG1	IRF4::PLAG1	chr14	23289056	23289057	rs17122779	GM128781:3010440	IRF4::PLAG1	75	44	chr14	23288990	23289065	23288960	23289094	IRF4	23288990	23289004	PLAG1	23289048	23289065	GM12878	loop_loop	2	IRF4::PLAG1+GM12878+2	0	TRUE	A_scrB_seq	TTGACTCTGAGGAGCAGAAAAAGCAGAGGGTAACGGGGACTGAAATGGAAAACAGGAAAACACATGTGTGACTCTGGCCGGCAGCGGAGGCCGGAGGACGCGGGGGGTGGAGCTGGGACTGCCTGGTCAGACCCA	IRF4::PLAG1+GM12878+2+A_scrBchr14	23289048	23289065	PLAG1	IRF4::PLAG1	chr14	23289056	23289057	rs17122779	GM12878and the second one PLAG1 was the one with the SNP rs17122779 (GM12878)<jupyter_code>mpra_counts = pd.read_csv('../../miseq/mpraanalyze_results/dk_norm/count_matrix_df_normlog.csv')
mpra_counts
km_df = mpra_counts[mpra_counts['Unnamed: 0'].str.contains('KLF5::MAFG\+GDSD6\+2')].melt('Unnamed: 0')
km_df.columns = ['id','sample','value']
km_df[['vocab_pair','tissue','instance','orient_bar']]= km_df['id'].str.split('+',expand=True)
km_df[['orient','barcode']]= km_df['orient_bar'].str.split('.',expand=True)
km_df[['mpra_tissue','rep']]= km_df['sample'].str.split('_',expand=True)
km_df.sort_values('orient',inplace=True)
km_df
ax = sns.boxplot(x="mpra_tissue", y="value", hue="orient",hue_order = ['scrA_scrB','scrA_B','A_scrB','A_B'],

                   data=km_df, palette="Set2", dodge=True)
ax.set_title('KLF5::MAFG+GDSD6+2')
ax = sns.boxplot(x="orient", y="value",# hue="orient",hue_order = ['scrA_scrB','scrA_B','A_scrB','A_B'],

                   data=km_df[km_df.mpra_tissue=='KCD6'], palette="Set2", dodge=True)
sns.swarmplot(x="orient", y="value",ax=ax,# hue="orient",hue_order = ['scrA_scrB','scrA_B','A_scrB','A_B'],
                   data=km_df[km_df.mpra_tissue=='GM'], color='black',dodge=True)
ax.set_title('KLF5::MAFG+GDSD6+2')
plt.savefig('../data/processed/snp_motifs/KLF5::MAFG+GDSD6+2.pdf')
ip_df = mpra_counts[mpra_counts['Unnamed: 0'].str.contains('IRF4::PLAG1\+GM12878\+2')].melt('Unnamed: 0')
ip_df.columns = ['id','sample','value']
ip_df[['vocab_pair','tissue','instance','orient_bar']]= ip_df['id'].str.split('+',expand=True)
ip_df[['orient','barcode']]= ip_df['orient_bar'].str.split('.',expand=True)
ip_df[['mpra_tissue','rep']]= ip_df['sample'].str.split('_',expand=True)
ip_df.sort_values('orient',inplace=True)
ip_df
ax = sns.boxplot(x="mpra_tissue", y="value", hue="orient",hue_order = ['scrA_scrB','scrA_B','A_scrB','A_B'],
                   data=ip_df, palette="Set2", dodge=True)
ax.set_title('IRF4::PLAG1+GM12878+2')
ax = sns.boxplot(x="orient", y="value",# hue="orient",hue_order = ['scrA_scrB','scrA_B','A_scrB','A_B'],

                   data=ip_df[ip_df.mpra_tissue=='GM'], palette="Set2", dodge=True)
sns.swarmplot(x="orient", y="value",ax=ax,# hue="orient",hue_order = ['scrA_scrB','scrA_B','A_scrB','A_B'],
                   data=ip_df[ip_df.mpra_tissue=='GM'], color='black',dodge=True)
ax.set_title('IRF4::PLAG1+GM12878+2')
plt.savefig('../data/processed/snp_motifs/IRF4::PLAG1+GM12878+2.pdf')<jupyter_output><empty_output><jupyter_text>Try 2: snps from cancer stuff but includes all things (28 hits)<jupyter_code>snp_df = pd.read_csv('../data/external/GWAS/snp_list_combined.csv',index_col=0)
snp_df.columns = ['Linked_SNP','Chr','Position','disease']
snp_bed_df = snp_df 
snp_bed_df['stop'] = snp_bed_df.Position+1
snp_bed_df = snp_bed_df[['Chr','Position','stop','Linked_SNP']]
snp_bed_df.columns = ['chr','start','stop','name']



snp_bed = pybedtools.BedTool.from_dataframe(snp_bed_df)
print(snp_df.shape)
display(snp_df[:5])


snp_vocab_all = get_snp_vocab(snp_bed)

snp_vocab_all = snp_vocab_all.merge(snp_df[['Linked_SNP','disease']]).drop_duplicates()
print(snp_vocab_all.shape)
snp_vocab_all<jupyter_output>(35203, 5)<jupyter_text>Try 3: noncancer linked snps<jupyter_code>snp_loc_df = pd.read_csv('../data/external/GWAS/causalSNP-combined-location.vcf',header=None, sep='\t')
snp_loc_df.columns = ['Chr','Position','Linked_SNP','Ref','Alt']
snp_loc_df['Chr'] = 'chr'+snp_loc_df.Chr.map(str)

snp_df = pd.read_csv('../data/external/GWAS/causalSNP-leadSNP-nonCancer-Oct23.txt',header=None, sep='\t')
snp_df.columns = ['Linked_SNP','beta','Index_SNP','pop','r2','Index_SNP2','disease']
snp_df = snp_df.merge(snp_loc_df, how='inner',on='Linked_SNP')
snp_bed_df = snp_df 
snp_bed_df['stop'] = snp_bed_df.Position+1
snp_bed_df = snp_bed_df[['Chr','Position','stop','Linked_SNP']].drop_duplicates()
snp_bed_df.columns = ['chr','start','stop','name']



# snp_bed = pybedtools.BedTool.from_dataframe(snp_bed_df)
print(snp_df.shape, snp_bed_df.shape)
display(snp_df[:5])


snp_vocab_all = get_snp_vocab(snp_bed).drop_duplicates()

# snp_vocab_all = snp_vocab_all.merge(snp_df[['Linked_SNP']]).drop_duplicates()
print(snp_vocab_all.shape)
snp_vocab_all



snp_vocab_all.Linked_SNP.unique().shape
snp_vocab_all.to_csv('../data/processed/snp_motifs/snp_vocab_all.csv')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1669052732044858, 'fraction_numerical': 0.05774278215223097, 'mean_word_length': 7.408224674022066, 'pattern_counts': {'":': 0, '<': 20, '<?xml version=': 0, '>': 22, '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': '212509', 'license': 'MIT', 'n_tokens_mistral': 4016, 'n_tokens_neox': 3642, 'n_words': 400, 'path': 'notebooks/16_snp_invocab.ipynb', 'repo_name': 'mguo123/pan_omics', 'size': '8382'}