Skip to content

Commit

Permalink
issues fixed
Browse files Browse the repository at this point in the history
  • Loading branch information
hajia018 committed Jun 23, 2023
1 parent 9d09f52 commit 7e2c174
Show file tree
Hide file tree
Showing 4 changed files with 18 additions and 15 deletions.
4 changes: 2 additions & 2 deletions corefuns/fdrsampleperm.py
Original file line number Diff line number Diff line change
Expand Up @@ -124,8 +124,8 @@ def fdrsampleperm(ssmFile, BPMindFile, pcut, minpath, N, genesetname=None):

# calling calculate_fdr() function to compute FDRs
fdrBPM1, fdrBPM2 = calculate_fdr(bpmdf, bpm_cols, bpm_pvdf, bpm_pv_cols, pcut, N, 'bpm')
fdrWPM1, fdrWPM2 = calculate_fdr(wpmdf, wpm_cols, wpm_pvdf, wpm_pv_cols, pcut, N, 'wpm')
fdrPATH1, fdrPATH2 = calculate_fdr(pathdf, path_cols, path_pvdf, path_pv_cols, pcut, N, 'path')
fdrWPM1, fdrWPM2 = calculate_fdr(wpmdf, wpm_cols, wpm_pvdf, wpm_pv_cols, 0.05, N, 'wpm')
fdrPATH1, fdrPATH2 = calculate_fdr(pathdf, path_cols, path_pvdf, path_pv_cols, 0.05, N, 'path')

bpm_ranksum = bpmdf[bpm_cols[0:1]]
bpm_ranksum.columns = ['bpm_ranksum']
Expand Down
9 changes: 5 additions & 4 deletions corefuns/pathway_map.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,14 +5,15 @@
import networkx as nx
import matplotlib.pyplot as plt
from matplotlib.lines import Line2D
from matplotlib import rcParams
from classes import fdrresultsclass
from corefuns import check_BPM_WPM_redundancy as cbwr
from networkx.drawing.nx_agraph import to_agraph
from matplotlib.backends.backend_pdf import PdfPages


def draw_map(project_dir,fdrcut,resultsfile,BPM_group_tmp,WPM_group_tmp,PATH_group_tmp,bpm_limit=20):

rcParams['font.family'] = 'sans-serif'
# load results file(FDRs)
#fdr_file = project_dir + '/' + resultsfile
pklin = open(resultsfile,'rb')
Expand Down Expand Up @@ -217,8 +218,8 @@ def draw_map(project_dir,fdrcut,resultsfile,BPM_group_tmp,WPM_group_tmp,PATH_gro
if val == 1:
p1 = used_pathways[i]
p2 = used_pathways[j]
G.add_edge(inds[p1],inds[p2],color='yellow')
G.add_edge(inds[p2],inds[p1],color='yellow')
G.add_edge(inds[p1],inds[p2],color='orange')
G.add_edge(inds[p2],inds[p1],color='orange')
elif val == -1:
p1 = used_pathways[i]
p2 = used_pathways[j]
Expand All @@ -244,7 +245,7 @@ def draw_map(project_dir,fdrcut,resultsfile,BPM_group_tmp,WPM_group_tmp,PATH_gro
legend_elements = [Line2D([0], [0], marker='o', color='w', label='pathway',markerfacecolor='lightblue', markersize=15)]
#legend_elements.append(Line2D([0], [0], marker='o', color='w', label='protective PATH pwathway',markerfacecolor='yellow', markersize=15))
#legend_elements.append(Line2D([0], [0], marker='o', color='w', label='risk PATH pwathway',markerfacecolor='blue', markersize=15))
legend_elements.append(Line2D([0], [0], color='yellow', lw=4, label='protective interaction'))
legend_elements.append(Line2D([0], [0], color='orange', lw=4, label='protective interaction'))
legend_elements.append(Line2D([0], [0], color='blue', lw=4, label='risk interaction'))
legend_elements.append(Line2D([0], [0], marker='o', color='w', label='self loops indicate WPM interactions',markerfacecolor='w', markersize=1))
plt.legend(handles=legend_elements, loc='lower left')
Expand Down
18 changes: 9 additions & 9 deletions scripts/data_checkpopulation.sh
Original file line number Diff line number Diff line change
Expand Up @@ -102,9 +102,9 @@ grep -w ASW ${popIDFile} | awk '{print $1 "\t" $2}' >> intermediate/sublist.tmp
grep -w YRI ${popIDFile} | awk '{print $1 "\t" $2}' >> intermediate/sublist.tmp

# add ${prj1000Pop} if it's not in CEU, CHB, YRI and ASW
if [[ ! "${prj1000Pop}" =~ ^(CEU|CHB|ASW|YRI)$ ]]; then
grep -w ${prj1000Pop} ${popIDFile} | awk '{print $1 "\t" $2}' >> intermediate/sublist.tmp
fi
#if [[ ! "${prj1000Pop}" =~ ^(CEU|CHB|ASW|YRI)$ ]]; then
# grep -w ${prj1000Pop} ${popIDFile} | awk '{print $1 "\t" $2}' >> intermediate/sublist.tmp
#fi

for pop in "${prj1000Pop[@]}"; do
if [[ ! "${pop}" =~ ^(CEU|CHB|ASW|YRI)$ ]]; then
Expand All @@ -122,13 +122,13 @@ plink --bfile intermediate/prj1000_tmp1 --keep intermediate/sublist.tmp --allow-
plink --bfile intermediate/prj1000_tmp2 --biallelic-only strict --allow-no-sex --make-bed --out intermediate/prj1000_tmp3 > /dev/null

# merge 1000 project data with study data
plink --bfile intermediate/prj1000_tmp3 --bmerge raw/${plinkFile} --allow-no-sex --make-bed --out intermediate/allpop_tmp > /dev/null
plink --bfile intermediate/prj1000_tmp3 --bmerge raw/${plinkFile} --allow-no-sex --make-bed --out intermediate/allpop_tmp > /dev/null 2> /dev/null

# if any SNPs in two datasets have different alleles, there will be error "variants with 3+ alleles present"
if [ -f "allpop_tmp-merge.missnp" ]; then
plink --bfile intermediate/prj1000_tmp3 --exclude intermediate/allpop_tmp-merge.missnp --allow-no-sex --make-bed --out intermediate/prj1000_tmp4 > /dev/null
plink --bfile raw/${plinkFile} --exclude intermediate/allpop_tmp-merge.missnp --allow-no-sex --make-bed --out intermediate/${plinkFile}_tmp1 > /dev/null
plink --bfile intermediate/prj1000_tmp4 --bmerge intermediate/${plinkFile}_tmp1 --allow-no-sex --make-bed --out intermediate/allpop_tmp > /dev/null
plink --bfile intermediate/prj1000_tmp4 --bmerge intermediate/${plinkFile}_tmp1 --allow-no-sex --make-bed --out intermediate/allpop_tmp > /dev/null 2> /dev/null
fi


Expand All @@ -149,21 +149,21 @@ z=0
#echo `echo CEU CHB ASW YRI ${prj1000Pop} | xargs -n1 | sort -u | xargs` StudyPop > intermediate/allpopid.txt

valid_values=("CEU" "CHB" "ASW" "YRI")
pop_values=("${valid_values[@]}" ${prj1000Pop[@]})
pop_values=("${valid_values[@]}" "${prj1000Pop[@]}")

for i in $(printf '%s\n' "${pop_values[*]}" | sort -u); do
for i in $(printf '%s\n' "${pop_values[*]}"); do
z+=" 0"
done

echo $(printf '%s\n' "${pop_values[*]}" | sort -u) StudyPop > intermediate/allpopid.txt
echo $(printf '%s\n' "${pop_values[@]}") StudyPop > intermediate/allpopid.txt


num_columns=$(( ${#pop_values[@]} + 1 ))
while read -r line; do
pattern=$(echo "$line" | awk '{print $1 " " $2}')
pop=$(grep -w "$pattern" "$popIDFile" | awk '{print $3}')

if [[ "${pop_values[@]}" =~ "$pop" ]]; then
if [[ -n "$pop" ]] && [[ "${pop_values[@]}" =~ "$pop" ]]; then
idx=$(printf '%s\n' "${pop_values[@]}" | awk -v p="$pop" 'BEGIN{ORS="\n";NR=0}$1==p{print NR}')
echo "$pattern $(echo "$z" | awk -v i="$idx" '{ $i=1; print }')"
else
Expand Down
2 changes: 2 additions & 0 deletions scripts/plotmds.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,7 @@
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt
from matplotlib import rcParams
import sys


Expand All @@ -27,6 +28,7 @@


if __name__ == '__main__':
rcParams['font.family'] = 'sans-serif'
mds_file = sys.argv[1]
pop_id_file = sys.argv[2]
prj1000Pop = sys.argv[3]
Expand Down

0 comments on commit 7e2c174

Please sign in to comment.