From 9a85ad752194846f26f2555d1f059f76ed31c43d Mon Sep 17 00:00:00 2001 From: Lan Hui Date: Wed, 29 Oct 2025 18:00:27 +0800 Subject: Also look for a run's library strategy/source in EXPERIMENT_LIBRARY_INFO_FILE --- Code/configure.py | 2 ++ Code/download_and_map.py | 22 +++++++++++++++++++--- Code/download_ena_experiment_records.py | 17 ++++++++--------- 3 files changed, 29 insertions(+), 12 deletions(-) diff --git a/Code/configure.py b/Code/configure.py index be04994..219ba1d 100644 --- a/Code/configure.py +++ b/Code/configure.py @@ -43,6 +43,8 @@ MIN_RNA_SEQ_INCREASE = -999 # minimum RNA-seq experiments needed when updating p UPDATE_NETWORK_LOG_FILE = '../Data/log/update.network.log.txt' # network update log. We should check this file from time to time. NEW_OR_UPDATED_CHIP_FILE = '../Data/log/new.or.updated.chip.file.txt' +EXPERIMENT_INFO_DIR = '../Data/information/EXPERIMENT_SET' +EXPERIMENT_LIBRARY_INFO_FILE = '../Data/information/experiment_library_info.txt' # is a run an RNA-Seq? Get the info from the xml file path EXPERIMENT_SET/EXPERIMENT/DESIGN/LIBRARY_DESCRIPTOR RNA_SEQ_INFO_DATABASE = '../Data/information/rnaseq_info_database.txt' # same as RNA_SEQ_INFO_FILE RNA_SEQ_INFO_DATABASE_JSON = '../Data/information/rnaseq_info_database.json' diff --git a/Code/download_and_map.py b/Code/download_and_map.py index e60b150..e9fb14b 100644 --- a/Code/download_and_map.py +++ b/Code/download_and_map.py @@ -23,7 +23,7 @@ import json from datetime import datetime ########################################################################################## -from configure import DAILY_MAP_NUMBER, MIN_FASTQ_FILE_SIZE, RNA_SEQ_INFO_FILE, DOWNLOADED_SRA_ID_LOG_FILE, IGNORED_SRA_ID_LOG_FILE, UPDATE_NETWORK_LOG_FILE, MAPPED_RDATA_DIR, RAW_RDATA_DIR, SALMON_MAP_RESULT_DIR +from configure import DAILY_MAP_NUMBER, MIN_FASTQ_FILE_SIZE, RNA_SEQ_INFO_FILE, EXPERIMENT_LIBRARY_INFO_FILE, DOWNLOADED_SRA_ID_LOG_FILE, IGNORED_SRA_ID_LOG_FILE, UPDATE_NETWORK_LOG_FILE, MAPPED_RDATA_DIR, RAW_RDATA_DIR, SALMON_MAP_RESULT_DIR ########################################################################################## def glob_files(directory, pattern): @@ -74,7 +74,7 @@ def get_list(fname): return result # only return unique elements -def make_download_list(mapped_dir, rna_data_info_dict): +def make_download_list(mapped_dir, rna_data_info_dict, experiment_library_info_dict={}): ''' Make next n sample IDs. These samples must have not been downloaded yet. @@ -89,6 +89,10 @@ def make_download_list(mapped_dir, rna_data_info_dict): small_ids = get_list(IGNORED_SRA_ID_LOG_FILE) # these data files are too small for run_id in sorted(rna_data_info_dict.keys(), reverse=True): # SRR first, then ERR, then DRR include_me_because_i_am_rnaseq = True if rna_data_info_dict[run_id]['library_strategy'].lower() == 'rna-seq' and rna_data_info_dict[run_id]['library_source'].lower() == 'transcriptomic' else False + associated_experiment_id = rna_data_info_dict[run_id]['experiment_id'] + if associated_experiment_id in experiment_library_info_dict: + is_rnaseq = 'rna-seq' in experiment_library_info_dict[associated_experiment_id]['library_strategy'].lower() and 'transcriptomic' in experiment_library_info_dict[associated_experiment_id]['library_source'] + include_me_because_i_am_rnaseq = include_me_because_i_am_rnaseq or is_rnaseq if not (run_id + '_quant.txt') in mapped_files and (not run_id in result) and (not run_id in small_ids) and (not run_id in mapped_run_ids) and include_me_because_i_am_rnaseq: # not mapped yet and is RNA-seq result.append(run_id) return result @@ -352,6 +356,17 @@ def read_ena_data_info_json(fname): return json_dict +def read_experiment_library_information(fname): + result = {} + with open(fname) as f: + for line in f: + line = line.strip() + lst = line.split('\t') + k = lst[0] + result[k] = {'library_strategy':lst[1], 'library_source':lst[2]} + return result + + def read_run_ids_from_file(fname): f = open(fname) lst = [] @@ -411,6 +426,7 @@ if not last_session_finished(DOWNLOADED_SRA_ID_LOG_FILE): # last session not fin sys.exit() rna_data_info_dict = read_ena_data_info_json(RNA_SEQ_INFO_FILE) # rna_data_info_dict contains only RNA-seq IDs. +experiment_library_info_dict = read_experiment_library_information(EXPERIMENT_LIBRARY_INFO_FILE) # Generate DRR/ERR/SRR ids to download if len(sys.argv) > 1: # user has provided a list of IDs in a file @@ -418,7 +434,7 @@ if len(sys.argv) > 1: # user has provided a list of IDs in a file DAILY_MAP_NUMBER = len(download_list) else: print('[download_and_map.py] Prepare download list ...') - download_list = make_download_list(MAPPED_RDATA_DIR, rna_data_info_dict) + download_list = make_download_list(MAPPED_RDATA_DIR, rna_data_info_dict, experiment_library_info_dict) print('[download_and_map.py] There are %d run IDs from which you could select %d of them.' % (len(download_list), DAILY_MAP_NUMBER)) # Make a record in log.txt diff --git a/Code/download_ena_experiment_records.py b/Code/download_ena_experiment_records.py index 14a4a96..c897e9c 100644 --- a/Code/download_ena_experiment_records.py +++ b/Code/download_ena_experiment_records.py @@ -18,6 +18,7 @@ import urllib.request import os, time, glob from parse_ena_xml import parse_experiment +from configure import EXPERIMENT_INFO_DIR def get_xml_for_experiment_id(eid): url = 'https://www.ebi.ac.uk/ena/browser/api/xml/%s' % (eid) @@ -26,23 +27,21 @@ def get_xml_for_experiment_id(eid): content = response.read() return content.decode('utf-8') - if __name__ == '__main__': - experiment_info_dir = '../Data/information/EXPERIMENT_SET' - if not os.path.exists(experiment_info_dir): - os.mkdir(experiment_info_dir) - already_downloaded_lst = [os.path.split(path)[1] for path in sorted(glob.glob(os.path.join(experiment_info_dir, '*')))] + if not os.path.exists(EXPERIMENT_INFO_DIR): + os.mkdir(EXPERIMENT_INFO_DIR) + already_downloaded_lst = [os.path.split(path)[1] for path in sorted(glob.glob(os.path.join(EXPERIMENT_INFO_DIR, '*')))] print('%d experiment info has already been downloaded' % (len(already_downloaded_lst))) - with open('../Data/information/experiment_ids_lacking_strategy_or_source.txt') as f: # the file experiment_ids_lacking_strategy_or_source.txt is produced by parse_ena_records.py + with open('../Data/information/experiment_ids_lacking_strategy_or_source.txt') as f: # the file experiment_ids_lacking_strategy_or_source.txt is produced by parse_ena_xml.py for line in f: experiment_id = line.strip() if experiment_id not in already_downloaded_lst: print(f'Downloading EXPERIMENT info for {experiment_id}') xml_content = get_xml_for_experiment_id(experiment_id) - fname = os.path.join(experiment_info_dir, experiment_id) + fname = os.path.join(EXPERIMENT_INFO_DIR, experiment_id) with open(fname, 'w', encoding='utf-8') as f2: f2.write(xml_content) d = parse_experiment(fname) - print(d[experiment_id]['library_strategy']) - print(d[experiment_id]['library_source']) + print(' ' + d[experiment_id]['library_strategy']) + print(' ' + d[experiment_id]['library_source']) time.sleep(6) -- cgit v1.2.1