-
Notifications
You must be signed in to change notification settings - Fork 1
Expand file tree
/
Copy pathmake_GOplots_dataframes.py
More file actions
242 lines (216 loc) · 12.5 KB
/
Copy pathmake_GOplots_dataframes.py
File metadata and controls
242 lines (216 loc) · 12.5 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
##############################################################################
### ###
### AUTHOR: CONTESSA A. RICCI, PhD ###
### MANUSCRIPT: DYSREGULATION OF MITOCHONDRIA-MEDIATED MATERNAL-FETAL ###
### INTERACTIONS IN HYPERTENSIVE DISORDERS OF PREGNANCY ###
### DOI: (TBD) ###
### STUDY PURPOSE: Reanalysis of longitudinal maternal RNAseq data from ###
### peripheral blood plasma and endpoint fetal RNAseq data ###
### from placenta at delivery. Goal is to understand ###
### consequences of mitochondrial gene dysregulation by ###
### examining expression patterns of genes the dysregulated ###
### mitochondrial genes are known to interact with ###
### DATA: GEO accession GSE154377 (maternal longitudinal data) ###
### GEO accession GSE114691 (fetal placental) ###
### Accessible via NCBI GEO datasets ###
### ###
### NOTES: All scripts assume all data has been compiled and is in correct ###
### format. This script utilizes some files that have been compiled ###
### either manually, in Python, or in R (or a combination thereof). ###
### ###
##############################################################################
# PURPOSE OF SCRIPT: Create dataframes for GOplot plots in R. Z-scores will be calculated in R.
#
#
# notes:
# - will use mtDEG-INTXs expression files generated by find_mtDEG-INTXs.py script:
# - maternal_del_mtDEG-INTXs_expression.csv
# - fetal_del_mtDEG-INTXs_expression.csv
# - will use GO enrichment files generated using GO enrichment analysis feature at geneontology.org:
# - maternal_delivery_BPresults.txt (analysis conducted on 12/06/2021)
# - maternal_delivery_CCresults.txt (analysis conducted on 12/06/2021)
# - fetal_delivery_BPresults.txt (analysis conducted on 12/06/2021)
# - fetal_delivery_CCresults.txt.txt (analysis conducted on 12/06/2021)
# note: Maternal and fetal GO enrichments were determined in separate analyses. GO enrichment
# categories (biological process, cellular component) were also determined in separate analyses.
# GO enrichment files were generated as follows:
# - Uniprot IDs outputted using find_mtDEG-INTXs.py script
# (fetal_del_mtDEG-INTXs_UniprotIDs.txt, maternal_del_mtDEG-INTXs_UniprotIDs.txt)
# were pasted into "Go Enrichment Analysis" dialog box.
# - For biological process enrichments "biological process" was selected from the
# dropdown menu immediately below dialog box where Uniprot IDs were pasted.
# Likewise, "cellular component" was selected for cellular component enrichments.
# - "Homo sapiens" was was specified as species using dropdown menu immediately below
# dropdown menu for GO enrichment categories.
# - Click the "Launch" button (located on same row as dropdown menu for species
# specification) to run analysis
# - Results for biological process or cellular component analysis were exported as a
# tab-delimited file by clicking the blue "Table" button, which is the first export
# option above the results (maternal_delivery_BPresults.txt,
# maternal_delivery_CCresults.txt, fetal_delivery_BPresults.txt,
# fetal_delivery_CCresults.txt)
# - will use Amigo database (http://amigo.geneontology.org/amigo):
# - amigo_db.txt (downloaded on 05/12/2021)
# note: will be used to match GO terms with maternal and fetal mtDEG-INTXs expression. Amigo file has entire
# entire Amgio database, comprised of GO accession IDs and terms for biological process GO terms, cellular
# component GO terms, and molecular function GO terms. Amigo file was generated as follows:
# - Create blank text file (amigo_db.txt)
# - Navigate to drill-down browser (http://amigo.geneontology.org/amigo/dd_browse)
# - Click on grey bubble to right of GO term category ("biological_process",
# "cellular_component", "molecular_function") with total terms in white text.
# - Click "retrieve gene products annotated to this term for this filter set"
# (hyperlinked text) at bottom of dialog box that pops up
# - Specify "Homo sapiens" under "Organism"
# - Click on "Custom DL (up to 1000000)" button located above terms listed
# - Move default options under "Selected fields" to "Available pool". Move
# "Acc (id)" and "Direct annotation (annotation_class_list)" from "Available pool"
# to "Selected fields"
# - Click "Download" at bottom of dialog box
# - Copy and paste all terms to amigo_db.txt text file
#
# OVERVIEW OF PROCESS:
# 1) COLLECT ALL UNIPROT IDs ASSOCIATED WITH ENRICHED GO TERMS
# 2) MATCH ENRICHED GO TERMS TO mtDEG-INTXs IDENTIFIED
# 3) MAKE PLOTTING FILES (MATERNAL AND FETAL AT DELIVERY)
###################################################
### ###### READ IN FILES ###### ###
###################################################
## READ IN AMIGO DATABASE ##
# [0] = Uniprot ID (in format UniProtKB:XXX)
# [1] = Gene name
# [2] = GO accession IDs associated with Uniprot ID, pipe ("|") delimited
amigo = open("amigo_db.txt", "r").read().split("\n")
## READ IN mtDEG-INTXs EXPRESSION FILES ##
# [0] = Ensembl transcript ID
# [1] = baseMean
# [2] = log2FoldChange
# [3] = lfcSE
# [4] = stat
# [5] = pvalue
# [6] = padj
# [7] = Uniprot_ID
maternal_expression = open("maternal_del_mtDEG-INTXs_expression.csv", "r").read().split("\n")
# [0] = Gene name
# [1] = DESeq2 (BH adjusted p-value)
# [2] = edgeR (BH adjusted p-value)
# [3] = ALDEx2 (effect size)
# [4] = Direction of expression
# [5] = Uniprot ID
fetal_expression = open("fetal_del_mtDEG-INTXs_expression.csv", "r").read().split("\n")
## READ IN GO ENRICHMENT FILES ##
# note: remove blank rows and informational rows
# [0] = GO term (term (GO:XXX))
# [1] = size of GO term (total genes belonging to GO term)
# [2] = total mtDEG-INTXs provided that belong to GO term
# [3] = proportion of provided mtDEG-INTXs that make up size of GO term
# [4] = fold enrichment
# [5] = raw P-value
# [6] = FDR (adjusted P-value)
maternal_BP = open("maternal_delivery_BPresults.txt", "r").read().split("\n")[11:]
maternal_CC = open("maternal_delivery_CCresults.txt", "r").read().split("\n")[11:]
fetal_BP = open("fetal_delivery_BPresults.txt", "r").read().split("\n")[11:]
fetal_CC = open("fetal_delivery_CCresults.txt", "r").read().split("\n")[11:]
################################################################
### ###### GO ENRICHMENTS DICTIONARIES ###### ###
################################################################
# make maternal GO enrichments dictionary
maternal_GOterms = {"BP":{}, "CC":{}}
for line in maternal_BP[1:]:
if len(line) > 1:
maternal_GOterms["BP"][line.split("\t")[0].split("(")[1][:-1]] = {"count":line.split("\t")[2],
"pval":line.split("\t")[6],
"GO_term":line.split("\t")[0].split("(")[0][:-1],
"Uniprot_IDs":[]}
for line in maternal_CC[1:]:
if len(line) > 1:
maternal_GOterms["CC"][line.split("\t")[0].split("(")[1][:-1]] = {"count":line.split("\t")[2],
"pval":line.split("\t")[6],
"GO_term":line.split("\t")[0].split("(")[0][:-1],
"Uniprot_IDs":[]}
# make fetal GO enrichments dictionary
fetal_GOterms = {"BP":{}, "CC":{}}
for line in fetal_BP[1:]:
if len(line) > 1:
fetal_GOterms["BP"][line.split("\t")[0].split("(")[1][:-1]] = {"count":line.split("\t")[2],
"pval":line.split("\t")[6],
"GO_term":line.split("\t")[0].split("(")[0][:-1],
"Uniprot_IDs":[]}
for line in fetal_CC[1:]:
if len(line) > 1:
fetal_GOterms["CC"][line.split("\t")[0].split("(")[1][:-1]] = {"count":line.split("\t")[2],
"pval":line.split("\t")[6],
"GO_term":line.split("\t")[0].split("(")[0][:-1],
"Uniprot_IDs":[]}
# populate GO enrichments dictionary with all possible Uniprot IDs associated with GO accession
for key, value in maternal_GOterms["BP"].items():
for line in amigo:
if key in line and "UniProtKB" in line.split("\t")[0]:
if line.split("\t")[0].split(":")[1] not in value["Uniprot_IDs"]:
value["Uniprot_IDs"].append(line.split("\t")[0].split(":")[1])
print("Maternal BP Uniprot IDs collected")
for key, value in maternal_GOterms["CC"].items():
for line in amigo:
if key in line and "UniProtKB" in line.split("\t")[0]:
if line.split("\t")[0].split(":")[1] not in value["Uniprot_IDs"]:
value["Uniprot_IDs"].append(line.split("\t")[0].split(":")[1])
print("Maternal CC Uniprot IDs collected")
for key, value in fetal_GOterms["BP"].items():
for line in amigo:
if key in line and "UniProtKB" in line.split("\t")[0]:
if line.split("\t")[0].split(":")[1] not in value["Uniprot_IDs"]:
value["Uniprot_IDs"].append(line.split("\t")[0].split(":")[1])
print("Fetal BP Uniprot IDs collected")
for key, value in fetal_GOterms["CC"].items():
for line in amigo:
if key in line and "UniProtKB" in line.split("\t")[0]:
if line.split("\t")[0].split(":")[1] not in value["Uniprot_IDs"]:
value["Uniprot_IDs"].append(line.split("\t")[0].split(":")[1])
print("Fetal CC Uniprot IDs collected")
########################################################
### ###### FINAL DICTIONARIES ###### ###
########################################################
# assign GO enrichment results to maternal mtDEG-INTXs identified
maternal_plotting_data = {"BP":{}, "CC":{}}
for line in maternal_expression[1:]:
for key, value in maternal_GOterms["BP"].items():
if len(line) > 1 and line.split(",")[7] in value["Uniprot_IDs"]:
maternal_plotting_data["BP"][line.split(",")[0] + "," + key] = {"Expression":line.split(",")[2],
"pval":value["pval"], "count":value["count"], "GO_term":value["GO_term"], "Uniprot_ID":line.split(",")[7]}
for key, value in maternal_GOterms["CC"].items():
if len(line) > 1 and line.split(",")[7] in value["Uniprot_IDs"]:
maternal_plotting_data["CC"][line.split(",")[0] + "," + key] = {"Expression":line.split(",")[2],
"pval":value["pval"], "count":value["count"], "GO_term":value["GO_term"], "Uniprot_ID":line.split(",")[7]}
# assign GO enrichment results to maternal mtDEG-INTXs identified
fetal_plotting_data = {"BP":{}, "CC":{}}
for line in fetal_expression[1:]:
for key, value in fetal_GOterms["BP"].items():
if len(line) > 1 and line.split(",")[5] in value["Uniprot_IDs"]:
fetal_plotting_data["BP"][line.split(",")[0] + "," + key] = {"Expression":line.split(",")[4],
"pval":value["pval"], "count":value["count"], "GO_term":value["GO_term"], "Uniprot_ID":line.split(",")[5]}
for key, value in fetal_GOterms["CC"].items():
if len(line) > 1 and line.split(",")[5] in value["Uniprot_IDs"]:
fetal_plotting_data["CC"][line.split(",")[0] + "," + key] = {"Expression":line.split(",")[4],
"pval":value["pval"], "count":value["count"], "GO_term":value["GO_term"], "Uniprot_ID":line.split(",")[5]}
########################################################
### ###### WRITE PLOTTING FILES ###### ###
########################################################
# note: z-scores are calculated in R. Column names meant to be compatible with GOplot package (R)
# maternal file
maternal_plotting_file = open("maternal_delivery_GOplot_data.csv", "w")
maternal_plotting_file.write("category,Ensembl ID,ID,term,count,genes,logFC,adj_pval\n")
for key, value in maternal_plotting_data["BP"].items():
maternal_plotting_file.write("BP" + "," + key + "," + value["GO_term"] +
"," + value["count"] + "," + value["Uniprot_ID"] + "," + value["Expression"] + "," + value["pval"] + "\n")
for key, value in maternal_plotting_data["CC"].items():
maternal_plotting_file.write("CC" + "," + key + "," + value["GO_term"] +
"," + value["count"] + "," + value["Uniprot_ID"] + "," + value["Expression"] + "," + value["pval"] + "\n")
# fetal file
fetal_plotting_file = open("fetal_delivery_GOplot_data.csv", "w")
fetal_plotting_file.write("category,Gene name,ID,term,count,genes,logFC,adj_pval\n")
for key, value in fetal_plotting_data["BP"].items():
fetal_plotting_file.write("BP" + "," + key + "," + value["GO_term"] +
"," + value["count"] + "," + value["Uniprot_ID"] + "," + value["Expression"] + "," + value["pval"] + "\n")
for key, value in fetal_plotting_data["CC"].items():
fetal_plotting_file.write("CC" + "," + key + "," + value["GO_term"] +
"," + value["count"] + "," + value["Uniprot_ID"] + "," + value["Expression"] + "," + value["pval"] + "\n")
print("Files written. Script complete")