### Install methylprep
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/general_walkthrough.html
Run this command in your terminal to install the methylprep package.
```bash
>>> pip install methylprep
```
--------------------------------
### Sample Sheet Metadata Info
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/general_walkthrough.html
Example output showing additional fields parsed from a sample sheet during the download process.
```text
INFO:methylprep.processing.pipeline:Found 5 additional fields in sample_sheet: source | histological diagnosis --> histological_diagnosis | gender | description | Sample_ID
```
--------------------------------
### View verbose CLI output for metadata download
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/faq.html
Example of the console output generated when running the metadata extraction command.
```text
INFO:methylprep.download.miniml:Downloading GSE84727_family.xml
INFO:methylprep.download.miniml:Downloaded GSE84727_family.xml
INFO:methylprep.download.miniml:Unpacking GSE84727_family.xml
INFO:methylprep.download.miniml:Downloaded and unpacked GSE84727
INFO:methylprep.download.miniml:MINiML file does not provide `methylprep_name` (sentrix_id_R00C00) for 847/847 samples.
INFO:methylprep.download.miniml:dropped Sentrix_ID (empty)
INFO:methylprep.download.miniml:dropped Sentrix_Position (empty)
INFO:methylprep.download.miniml:source == description; dropping source
INFO:methylprep.download.miniml:dropped `platform` ({'GPL13534'})
INFO:methylprep.download.miniml:title == Sample_Name; dropping title
INFO:methylprep.download.miniml:Final samplesheet contains 847 rows and 7 columns
INFO:methylprep.download.miniml:['GSM_ID', 'Sample_Name', 'sentrixids', 'Sex', 'age', 'disease_status', 'description']
```
--------------------------------
### Process multiple samples via CLI
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/cli.html
Example command for processing specific samples when a sample sheet is not provided.
```bash
python -m methylprep process -d . --all --no_sample_sheet -n Sample_1 Sample_2 Sample_3
```
--------------------------------
### Install methylsuite package
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/index.html
Use pip to install the complete methylsuite, which includes methylprep, methylcheck, and methylize.
```bash
>>> pip install methylsuite
```
--------------------------------
### Seek to Manifest Start
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/files/manifests.html
Positions the file pointer at the beginning of the data section by searching for the 'IlmnID' header.
```python
@staticmethod
def seek_to_start(manifest_file):
""" find the start of the data part of the manifest. first left-most column must be "IlmnID" to be found."""
reset_file(manifest_file)
current_pos = manifest_file.tell()
header_line = manifest_file.readline()
while not header_line.startswith(b'IlmnID'):
current_pos = manifest_file.tell()
if not header_line: #EOF
raise EOFError("The first (left-most) column in your manifest must contain 'IlmnID'. This defines the header row.")
header_line = manifest_file.readline()
if current_pos == 0:
manifest_file.seek(current_pos)
else:
manifest_file.seek(current_pos - 1)
```
--------------------------------
### Initialize directory structure
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/download/process_data.html
Creates necessary directories for dictionaries and beta values for each platform.
```python
def initialize(path):
"""Creates one directory for dictionaries and one directory for beta_values per platform
Arguments:
path [required]
the path to the directory to create the platform directories"""
if not Path(path).is_dir():
#LOGGER.debug(f"Created {path} directory.")
Path(path).mkdir(parents=True, exist_ok=True)
for platform in PLATFORMS:
#if not os.path.exists(f"{path}/{platform}_beta_values"):
# #LOGGER.debug(f"Created {platform} beta_values directory")
# os.mkdir(f"{path}/{platform}_beta_values")
if not os.path.exists(f"{path}/{platform}_dictionaries"):
#LOGGER.debug(f"Created {platform} dictionaries directory")
os.mkdir(f"{path}/{platform}_dictionaries")
```
--------------------------------
### initialize(path)
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/download/process_data.html
Creates one directory for dictionaries and one directory for beta_values per platform.
```APIDOC
## initialize
### Description
Creates one directory for dictionaries and one directory for beta_values per platform.
### Arguments
- **path** (required) - The path to the directory to create the platform directories.
```
--------------------------------
### Get export file path
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/models/samples.html
Determines the target path for exporting data as a CSV file.
```python
def get_export_filepath(self, extension='csv'):
""" Called by run_pipeline to find the folder/filename to export data as CSV, but CSV file doesn't exist yet."""
return self.get_filepath(extension, 'processed', verify=False)
```
--------------------------------
### Accessing methylprep CLI help
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/cli.html
Displays the main help menu for the methylprep package.
```bash
>>> python -m methylprep
usage: methylprep [-h] [-v] {process,sample_sheet} ...
Utility to process methylation data from Illumina IDAT files
positional arguments:
{process,sample_sheet}
process process help
sample_sheet sample sheet help
optional arguments:
-h, --help show this help message and exit
-v, --verbose Enable verbose logging
```
--------------------------------
### SNP Probe Notebook Logic
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/models/probes.html
Example logic for filtering SNP probes from a manifest dataframe.
```python
#manifest = data_containers[0].manifest.data_frame
#snpProbesI = manifest[manifest['probe_type']=='SnpI']
#snpProbesI_Grn = snpProbesI[snpProbesI['Color_Channel']=='Grn']
#snpProbesI_Red = snpProbesI[snpProbesI['Color_Channel']=='Red']
#snpProbesII = manifest[manifest['probe_type']=='SnpII']
```
--------------------------------
### Initialize Sample Sheet and Validate Samples
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/processing/pipeline.html
Handles sample sheet loading, field renaming logging, and validation of user-provided sample names against the sample sheet.
```python
if make_sample_sheet:
create_sample_sheet(data_dir)
try:
sample_sheet = get_sample_sheet(data_dir, filepath=sample_sheet_filepath)
except Exception as e:
# e will be 'Too many sample sheets in this directory.'
instructions = check_array_folders(data_dir, verbose=True) # prints instructions for GEO multi-array data packages.
if instructions != []:
instructions = '\n'.join(instructions)
print(f"This folder contains idats for multiple types of arrays. Run each array separately:\n{instructions}")
sys.exit(0)
raise Exception(e)
samples = sample_sheet.get_samples()
if sample_sheet.renamed_fields != {}:
show_fields = []
for k,v in sample_sheet.renamed_fields.items():
if v != k:
show_fields.append(f"{k} --> {v}")
else:
show_fields.append(f"{k}")
LOGGER.info(f"Found {len(show_fields)} additional fields in sample_sheet:\n{' | '.join(show_fields)}")
if sample_name is not None:
if not isinstance(sample_name,(list,tuple)):
raise SystemExit(f"sample_name must be a list of sample_names")
matched_samples = [sample.name for sample in samples if sample.name in sample_name]
if set(matched_samples) != set(sample_name):
possible_sample_names = [sample.name for sample in samples]
unmatched_samples = [_sample for _sample in sample_name if _sample not in possible_sample_names]
raise SystemExit(f"Your sample_name filter does not match the samplesheet; these samples were not found: {unmatched_samples}")
```
--------------------------------
### Initialize a SampleSheet instance
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/files/sample_sheets.html
Generates a SampleSheet instance for a given directory. If no filepath is provided, it attempts to locate one automatically.
```python
def get_sample_sheet(dir_path, filepath=None):
"""Generates a SampleSheet instance for a given directory of processed data.
Arguments:
dir_path {string or path-like} -- Base directory of the sample sheet and associated IDAT files.
Keyword Arguments:
filepath {string or path-like} -- path of the sample sheet file if provided, otherwise
one will try to be found. (default: {None})
Returns:
[SampleSheet] -- A SampleSheet instance.
"""
LOGGER.debug('Reading sample sheet')
if not filepath:
filepath = find_sample_sheet(dir_path)
data_dir = PurePath(filepath).parent
return SampleSheet(filepath, data_dir)
```
--------------------------------
### Example Probe Sequence
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/introduction/introduction.html
A sample oligonucleotide probe sequence used to evaluate the methylation state of a CpG site.
```text
CTACAAATACGACACCCGCAACCCATATTTCATATATTATCTCATTTAAC
```
--------------------------------
### Get Data Types
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/files/manifests.html
Returns a dictionary mapping columns to their data types, specifically setting Address IDs to Int64.
```python
def get_data_types(self):
data_types = {
key: str for key in self.columns
}
data_types['AddressA_ID'] = 'Int64' #'float64' -- dtype found only in pandas 0.24 or greater
data_types['AddressB_ID'] = 'Int64' #'float64'
return data_types
```
--------------------------------
### Initialize and run the processing pipeline
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/processing/pipeline.html
Configures the processing steps, export options, and estimator for the pipeline execution.
```python
def make_pipeline(data_dir='.', steps=None, exports=None, estimator='beta', **kwargs):
"""Specify a list of processing steps for run_pipeline, then instantiate and run that pipeline.
steps:
list of processing steps
['all', 'infer_channel_switch', 'poobah', 'quality_mask', 'noob', 'dye_bias']
exports:
list of files to be saved; anything not specified is not saved; ['all'] saves everything.
['all', 'csv', 'poobah', 'meth', 'unmeth', 'noob_meth', 'noob_unmeth', 'sample_sheet_meta_data', 'mouse', 'control']
estimator:
which final format?
[beta | m_value | copy_number | None (returns containers instead)]"""
```
--------------------------------
### Example of numeric encoding in column headers
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/faq.html
This list represents a column header where 'control' is encoded as '1', which the parser cannot interpret.
```python
['GSM_ID', 'Sample_Name', 'diseasestatus (1=control, 2=scz patient)', 'source', 'gender', 'sample type', 'plate', 'sentrix barcode', 'sentrix position', 'well id', 'age', 'used_in_analysis', 'description']
```
--------------------------------
### Initialize Pipeline for GEO Datasets
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/download/geo.html
Sets up the configuration for running methylprep on GEO datasets, optimizing for download speeds.
```python
[docs]def pipeline_find_betas_any_source(**kwargs):
"""beta_bake: Sets up a script to run methylprep that saves directly to path or S3.
The slowest part of processing GEO datasets is downloading, so this handles that.
```
--------------------------------
### Accessing process command help
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/cli.html
Displays the help menu and required arguments for the process command.
```bash
>>> python -m methylprep process
[Error]:
the following arguments are required: -d/--data_dir
usage: methylprep process [-h] -d DATA_DIR
[--array_type {custom,27k,450k,epic,epic+,mouse}]
[-m MANIFEST] [-s SAMPLE_SHEET] [--no_sample_sheet]
[-n [SAMPLE_NAME [SAMPLE_NAME ...]]] [-b] [-v]
[--batch_size BATCH_SIZE] [-u] [-e] [-x]
[-i {float64,float32,float16}] [-c] [--poobah]
[--export_poobah] [--minfi] [--no_quality_mask] [-a]
Process Illumina IDAT files, producing NOOB, beta-value, or m_value corrected
scores per probe per sample
```
--------------------------------
### Get Probe Details
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/files/manifests.html
Retrieves probe information based on type and optional color channel, used for mapping probes to names.
```python
def get_probe_details(self, probe_type, channel=None):
"""used by infer_channel_switch. Given a probe type (I, II, SnpI, SnpII, Control) and a channel (Channel.RED | Channel.GREEN),
this will return info needed to map probes to their names (e.g. cg0031313 or rs00542420), which are NOT in the idat files."""
if not isinstance(probe_type, ProbeType):
raise Exception('probe_type is not a valid ProbeType')
if channel and not isinstance(channel, Channel):
raise Exception('channel not a valid Channel')
data_frame = self.data_frame
probe_type_mask = data_frame['probe_type'].values == probe_type.value
if not channel:
return data_frame[probe_type_mask]
channel_mask = data_frame['Color_Channel'].values == channel.value
return data_frame[probe_type_mask & channel_mask]
```
--------------------------------
### Create a sample sheet
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/files/sample_sheets.html
Initializes the creation of a samplesheet.csv file from IDAT files found in a GEO series directory.
```python
def create_sample_sheet(dir_path, matrix_file=False, output_file='samplesheet.csv',
sample_type='', sample_sub_type=''):
"""Creates a samplesheet.csv file from the .IDAT files of a GEO series directory
```
--------------------------------
### Identify XML metadata parsing error
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/faq.html
Example of an incorrectly tagged characteristic in a MINiML file that causes samplesheet processing errors.
```xml
3998567027_R02C02
```
--------------------------------
### Build Sample Objects
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/files/sample_sheets.html
Processes the internal DataFrame to instantiate Sample objects.
```python
def build_samples(self):
"""Builds Sample objects from the processed sample sheet rows.
Added to Sample as class_method: if the idat file is not in the same folder, (check if exists!) looks recursively for that filename and updates the data_dir for that Sample.
"""
self.__samples = []
#logging.info('Building samples')
for _index, row in self.__data_frame.iterrows():
if self.alt_headers:
sentrix_id = row['SentrixBarcode_A'].strip()
sentrix_position = row['SentrixPosition_A'].strip()
else:
sentrix_id = row['Sentrix_ID'].strip()
sentrix_position = row['Sentrix_Position'].strip()
if not (sentrix_id and sentrix_position):
continue
sample = Sample(
data_dir=self.data_dir, # this assumes the .idat files are in the same folder with the samplesheet.
sentrix_id=sentrix_id,
sentrix_position=sentrix_position,
**row,
)
if sample.renamed_fields != {}:
self.renamed_fields.update(sample.renamed_fields)
self.fields.update(sample.fields)
self.__samples.append(sample)
```
--------------------------------
### Bisulfite Conversion Example
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/introduction/introduction.html
Visual representation of how bisulfite sequencing converts unmethylated cytosines to thymines while leaving methylated cytosines unchanged.
```text
Original DNA sequence: A C G A T C G C G A
Bisulfite converted sequence: A C G A T T G T G A
```
--------------------------------
### Read SNP Probes
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/files/manifests.html
Reads SNP probes from a manifest file, filtering for those starting with 'rs' and ensuring address IDs are float64.
```python
def read_snp_probes(self, manifest_file):
""" Unlike cpg and control probes, these rs probes are NOT sequential in all arrays. """
#LOGGER.info(f'Reading snp probes: {Path(manifest_file.name).stem} --> {snp_df.shape[0]} found')
self.seek_to_start(manifest_file)
# since these are not sequential, loading everything and filtering by IlmnID.
snp_df = pd.read_csv(
manifest_file,
low_memory=False)
# 'O' type columns won't match in SigSet, so forcing float64 here. Also, float32 won't cover all probe IDs; must be float64.
snp_df = snp_df[snp_df['IlmnID'].str.match('rs', na=False)].astype({'AddressA_ID':'float64', 'AddressB_ID':'float64'})
return snp_df
```
--------------------------------
### run_series(id, path, dict_only=False, batch_size=BATCH_SIZE, clean=True, abort_if_no_idats=True, decompress=True)
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/download/process_data.html
Downloads IDATs and metadata for a series, then generates a metadata dictionary and beta value matrix for each platform in the series.
```APIDOC
## run_series(id, path, dict_only=False, batch_size=BATCH_SIZE, clean=True, abort_if_no_idats=True, decompress=True)
### Description
Downloads the IDATs and metadata for a series then generates one metadata dictionary and one beta value matrix for each platform in the series.
### Parameters
- **id** (string) - Required - The series ID (can be a GEO or ArrayExpress ID).
- **path** (string) - Required - The path to the directory to download the data to.
- **dict_only** (boolean) - Optional - If True, downloads idat files and meta data and creates data dictionaries for each platform, but does not process them further.
- **batch_size** (integer) - Optional - The batch_size to use when processing samples (number of samples run at a time).
- **clean** (boolean) - Optional - If True, removes intermediate processing files.
- **abort_if_no_idats** (boolean) - Optional - If True, aborts if no IDAT files are found.
- **decompress** (boolean) - Optional - If True, decompresses downloaded files.
```
--------------------------------
### Build composite dataset with composite
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/special_cases.html
Download and process multiple datasets listed in a text file into a single directory.
```bash
>>> python -m methylprep composite -l ~/tutorial/glioblastoma.txt -d ~/tutorial
```
--------------------------------
### Configure Working Directory and Batch Processing
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/download/geo.html
Sets up the temporary working directory based on EFS or local data directory configurations and iterates through GEO IDs to process data.
```python
#1: set working folder, and make sure /mnt/efs exists if used.
if kwargs.get('efs') and kwargs['efs'] is not None:
EFS = kwargs['efs']
mounted = False
try:
#efs_files = list(Path(EFS).glob('*'))
mounted = Path(EFS).exists()
LOGGER.info(f"{EFS} exists: {mounted}")
except Exception as e:
LOGGER.error(f"{EFS} error: {e}")
if not mounted and os.name != 'nt':
result = subprocess.run(["df", "-h"], stdout=subprocess.PIPE)
LOGGER.warning(f"EFS mount [df -h]: {result.stdout.decode('utf-8')}")
raise FileNotFoundError("This batch function has no {EFS} mounted. Cannot proceed.")
working = tempfile.TemporaryDirectory(dir=EFS)
elif not kwargs.get('efs') and kwargs.get('data_dir') != None:
if not Path(kwargs.get('data_dir')).exists():
Path(kwargs.get('data_dir')).mkdir(parents=True, exist_ok=True)
working = tempfile.TemporaryDirectory(dir=kwargs['data_dir'])
EFS = working.name
elif not kwargs.get('efs') and kwargs.get('data_dir') == None:
working = tempfile.TemporaryDirectory(dir='.')
EFS = working.name
if kwargs.get('data_dir') == None:
kwargs['data_dir'] = '.' # CLI seems to pass in None and the get() doesn't get it.
LOGGER.info(f"Starting batch GEO pipeline processor for {geo_ids}.")
#1: run EACH geo series through downloader/processor
each_geo_result = {}
for geo_id in geo_ids:
# this download processed data, or idats, or -tbl-1.txt files inside the _family.xml.tgz files.
# also downloads meta_data
# download_geo_processed() downloads nothing if idats exist.
result = download_geo_processed(geo_id, working, verbose=kwargs.get('verbose', False), use_headers=True)
LOGGER.info(f"DEBUG: download_geo_processed result: {result}")
extracted_files = [_file.name for _file in list(Path(working.name).rglob('*')) if not _file.is_dir()]
zipfile_paths = [_file for _file in list(Path(working.name).rglob('*')) if not _file.is_dir()]
# LOGGER.info(f"DEBUG: zipfile_paths 465: {zipfile_paths}")
#2: CHECK: Are all files in working.dir or a sub-folder? Must return path to each file off of workdir.
zipfile_names = []
for zipfile_path in zipfile_paths:
#LOGGER.info(f"DEBUG: checking {zipfile_path.name}")
working_name_last_part = Path(working.name).parts[-1]
_file = Path(working.name, zipfile_path.name)
if _file.exists() and _file.is_file(): # no sub-path used
zipfile_names.append(zipfile_path.name)
#LOGGER.info(f"DEBUG: found {zipfile_path.name}")
elif Path(zipfile_path).exists() and Path(zipfile_path).is_file():
for part in reversed(Path(zipfile_path).parts):
if path == working_name_last_part:
new_zipfile_name = f"{working_name_last_part}/{zipfile_path.name}"
zipfile_names.append(new_zipfile_name)
LOGGER.info(f"Found file {zipfile_path.name} was in some sub-folder, so returning as {new_zipfile_name}.")
else:
#raise FileNotFoundError(f"ERROR: Downloaded a file but cannot parse the path to it: {zipfile_path}")
return {"filenames": zipfile_names, "tempdir": working, "error": f"Downloaded a file but cannot parse the path to it: {zipfile_path}"}
# LOGGER.info(f"DEBUG zipfile_names 485: {zipfile_names}")
#3: memory check / debug
if os.name != 'nt' and kwargs.get('verbose') == True:
total_disk_space = subprocess.check_output(['du','-sh', EFS]).split()[0].decode('utf-8')
LOGGER.info(f"Tempfolder {EFS} contains {len(list([str(k) for k in Path(EFS).rglob('*')]))} files, {total_disk_space} total.")
LOGGER.info(f"DEBUG: {EFS} all files: {list([str(k) for k in Path(working.name).rglob('*')])}")
#4: LOGGING OUT -- result dict tells this function what it found.
if result['found_idats'] == False and result['processed_files'] == True and result['tbl_txt_files'] == False:
LOGGER.info(f"Found {len(extracted_files)} files for {geo_id}.")
if result['found_idats'] == False and result['processed_files'] == False and result['tbl_txt_files'] == False:
LOGGER.warning(f"No downloadable methylation data found for {geo_id}.")
continue
if result['tbl_txt_files'] == True:
LOGGER.info(f"Found -tbl-1.txt files with meta data for {geo_id}.")
if result['tbl_txt_files'] == True:
LOGGER.info(f"Found processed csv files for {geo_id}.")
```
--------------------------------
### Download and extract MINiML files
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/download/miniml.html
Downloads files from an FTP server, tracks progress with tqdm, and extracts the target XML file from a tar archive.
```python
local_files = []
for miniml_filename in file_parts:
LOGGER.info(f"Downloading {miniml_filename}")
miniml_ftp_path = f"miniml/{miniml_filename}"
miniml_local_path = f"{series_path}/{miniml_filename}"
miniml_file = open(miniml_local_path, 'wb')
try:
# Alternatively, as a hacky workaround, you might try download a file -- any file -- before sending your SIZE command. With request.UseBinary = true for that request, it should cause your client to send the "TYPE I" command to the FTP server. (And it won't matter if that download request fails; the TYPE command will still have been sent.)
filesize = ftp.size(miniml_ftp_path)
with tqdm(unit = 'b', unit_scale = True, leave = False, miniters = 1, desc = geo_id, total = filesize) as tqdm_instance:
def tqdm_callback(data):
tqdm_instance.update(len(data))
miniml_file.write(data)
ftp.retrbinary(f"RETR {miniml_ftp_path}", tqdm_callback)
except Exception as e:
if not str(e).startswith('550'):
LOGGER.info(e)
if verbose:
LOGGER.info('Failed to create a progress bar, but it is downloading...')
ftp.retrbinary(f"RETR {miniml_ftp_path}", miniml_file.write)
miniml_file.close()
if not Path(miniml_local_path).exists():
LOGGER.error(f"{miniml_local_path} missing")
min_tar = tarfile.open(miniml_local_path)
for file in min_tar.getnames(): # this only extracts the one .xml file we need.
#LOGGER.info(f"DEBUG miniml_filename {expected_miniml_filename} == tar_file {file} | {file == expected_miniml_filename}")
if file == expected_miniml_filename:
min_tar.extract(file, path=series_path)
if Path(miniml_local_path).exists() and remove_tgz:
Path(miniml_local_path).unlink()
if verbose:
LOGGER.info(f"Downloaded and unpacked {miniml_filename}")
for extracted_file in Path(series_path).glob('*'):
if verbose:
LOGGER.info(f"Extracted {extracted_file}")
if expected_miniml_filename == extracted_file.name:
miniml_local_path = f"{series_path}/{extracted_file.name}"
if verbose:
LOGGER.info(f"*Identified* meta_data: {extracted_file}")
local_files.append(str(extracted_file))
ftp.quit()
return {'meta_data': miniml_local_path, 'local_files': local_files}
```
--------------------------------
### run_series_list(list_file, path, dict_only=False, batch_size=BATCH_SIZE, **kwargs)
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/download/process_data.html
Downloads IDATs and metadata for a list of series, creating metadata dictionaries and dataframes of sample beta_values.
```APIDOC
## run_series_list
### Description
Downloads the IDATs and metadata for a list of series, creating metadata dictionaries and dataframes of sample beta_values.
### Arguments
- **list_file** (required) - The name of the file containing a list of GEO_IDS and/or Array Express IDs to download and process. This file must be located in the directory data is downloaded to.
- **path** (required) - The path to the directory to download the data to.
- **dict_only** (optional) - If True, only downloads data and creates dictionaries for each platform.
- **batch_size** (optional) - The batch_size to use when processing samples (number of samples run at a time).
```
--------------------------------
### Download Default Manifest
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/files/manifests.html
Downloads the required manifest file if it is not already present on the local filesystem.
```python
@staticmethod
def download_default(array_type, on_lambda=False):
"""Downloads the appropriate manifest file if one does not already exist.
Arguments:
array_type {ArrayType} -- The type of array to process.
Returns:
[PurePath] -- Path to the manifest file.
"""
dir_path = Path(MANIFEST_DIR_PATH).expanduser()
if on_lambda:
dir_path = Path(MANIFEST_DIR_PATH_LAMBDA).expanduser()
filename = ARRAY_TYPE_MANIFEST_FILENAMES[array_type]
filepath = Path(dir_path).joinpath(filename)
if Path.exists(filepath):
return filepath
LOGGER.info(f"Downloading manifest: {Path(filename).stem}")
src_url = urljoin(MANIFEST_REMOTE_PATH, filename)
download_file(filename, src_url, dir_path)
return filepath
```
--------------------------------
### Manage Platform Metadata and Sample Sheets
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/download/miniml.html
Handles platform-specific metadata storage, merges sample sheets, and triggers the creation of local sample sheet files.
```python
if platform in geo_platforms:
meta_dicts[platform][accession] = attributes_dir
samples_dict[platform][accession] = title
elif warned_once is False:
LOGGER.warning(f"{platform} is not among the supported platforms {geo_platforms}")
warned_once = True
if missing_methylprep_names > 0:
LOGGER.info( f"MINiML file does not provide (Sentrix_id_R00C00) for {missing_methylprep_names}/{len(samples)} samples." )
for platform in geo_platforms:
if verbose:
LOGGER.info(f"Platform {platform} --> {len(meta_dicts[platform])} samples")
if meta_dicts[platform]:
#Path(f"{data_dir}/{platform}").mkdir(parents=True, exist_ok=True)
meta_dicts[platform] = meta_from_idat_filenames(data_dir, meta_dicts[platform])
if merge == True:
# find samplesheet.csv with idat ID/Position and merge in here.
meta_dicts[platform] = merge_sample_sheets(data_dir, meta_dicts[platform])
sample_sheet_from_miniml(geo_id, data_dir, platform, samples_dict[platform], meta_dicts[platform],
save_df=True, extract_controls=extract_controls, require_keyword=require_keyword, verbose=verbose)
if sync_idats:
paths = list(Path(data_dir).rglob(f'{geo_id}_{platform}_samplesheet.csv'))
if len(paths) > 0:
remove_idats_not_in_samplesheet(paths[0], data_dir)
else:
LOGGER.warning(f"Could not locate file ({f'{geo_id}_{platform}_samplesheet.csv'}) in path of {data_dir}.")
cleanup(data_dir)
return result['local_files']
```
--------------------------------
### Generate sample sheet from MINiML
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/download/miniml.html
Constructs a sample sheet dictionary by combining sample titles and metadata, specifically handling the mapping of GSM IDs to Sentrix information.
```python
def sample_sheet_from_miniml(geo_id, series_path, platform, samp_dict, meta_dict, save_df=False, extract_controls=False, require_keyword=None, verbose=True):
"""Creates a sample_sheet for all samples of a particular platform for a given series
Arguments:
geo_id [required]
the GEO Accension for the desired series
series_path [required]
the directory containing the series data
platform [required]
the platform to generate a sample sheet for
samp_dict
the dictionary of {sample GSM_ID: title} for the given platform
meta_dict
the dictionary of {sample GSM_ID: meta_dict} for the given platform
Notes:
cleans up redundant columns, including Sample_Name == title and platform
"""
# id and position are NOT in the miniml file, must match with GSMID later
_dict = {'GSM_ID': [], 'Sample_Name': [], 'Sentrix_ID': [], 'Sentrix_Position': []}
for GSM_ID, sample_title in samp_dict.items():
_dict['GSM_ID'].append(GSM_ID)
_dict['Sample_Name'].append(sample_title)
for key, value in meta_dict[GSM_ID].items():
if key not in _dict:
_dict[key] = []
_dict[key].append(value)
```
--------------------------------
### methylprep.make_pipeline
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/source/modules.html
Specify a list of processing steps for run_pipeline, then instantiate and run that pipeline.
```APIDOC
## methylprep.make_pipeline(data_dir='.', steps=None, exports=None, estimator='beta', **kwargs)
### Description
Specify a list of processing steps for run_pipeline, then instantiate and run that pipeline. This function feeds a class that runs the run_pipeline function of transforms with a final estimator.
### Parameters
- **data_dir** (str) - Optional - The directory containing the data. Defaults to '.'
- **steps** (list) - Optional - List of processing steps. Options: ['all', 'infer_channel_switch', 'poobah', 'quality_mask', 'noob', 'dye_bias']
- **exports** (list) - Optional - List of files to be saved. Options: ['all', 'csv', 'poobah', 'meth', 'unmeth', 'noob_meth', 'noob_unmeth', 'sample_sheet_meta_data', 'mouse', 'control']
- **estimator** (str) - Optional - The final format. Options: ['beta', 'm_value', 'copy_number', None]. Defaults to 'beta'.
- **kwargs** (dict) - Optional - Additional configuration settings including 'sesame', 'batch_size', 'bit', 'low_memory', 'poobah_decimals', 'poobah_sig', 'verbose', and 'debug'.
```
--------------------------------
### Generate Sample Sheet from IDAT Files
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/download/miniml.html
Parses IDAT files in a directory to create a sample sheet CSV and optionally a metadata pickle file.
```python
series_dir = Path(f"{series_path}/{platform}")
idat_files = series_dir.glob('*Grn.idat')
_dict = {'GSM_ID': [], 'Sample_Name': [], 'Sentrix_ID': [], 'Sentrix_Position': []}
for idat in idat_files:
filename = str(idat).split("/")[-1]
if re.match('(GSM[0-9]+_[0-9a-zA-Z]+_R0[0-9].0[0-9].Grn.idat)', filename):
split_filename = filename.split("_")
_dict['GSM_ID'].append(split_filename[0])
_dict['Sentrix_ID'].append(split_filename[1])
_dict['Sentrix_Position'].append(split_filename[2])
else:
raise ValueError(f"{filename} has an unexpected naming format")
# confusing logic here, names_dir is unordered, but retrieving by key in correct order alleviates
for key in _dict['GSM_ID']:
_dict['Sample_Name'].append(platform_samples_dict[key])
df = pd.DataFrame(data=_dict)
df.to_csv(path_or_buf=(Path(f"{series_path}/{platform}", 'samplesheet.csv')),index=False)
if save_df:
df.to_pickle(Path(f"{series_path}/{platform}", 'sample_sheet_meta_data.pkl'))
```
--------------------------------
### Convert MINiML via CLI
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/source/modules.html
Command line interface usage for parsing GEO XML files and generating sample sheet metadata without requiring IDAT files.
```bash
python -m meta_data -i GSExxxxx -d
```
--------------------------------
### Create a meta_data.pkl file
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/special_cases.html
Generates a meta_data.pkl file from a GEO dataset identifier without performing a full download or process.
```bash
>>> python -m methylprep meta_data -i GSE147391 -d ~/tutorial
```
--------------------------------
### Generate Sample Sheet from IDAT Files
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/files/sample_sheets.html
Parses a directory for IDAT files to construct a sample sheet CSV. Requires a valid directory path and optionally supports Series Matrix file integration for sample naming.
```python
sample_dir = Path(dir_path)
if not sample_dir.is_dir():
raise FileNotFoundError(f'{dir_path} is not a valid directory path')
idat_files = sample_dir.rglob('*Grn.idat*') #.gz OK
_dict = {'GSM_ID': [], 'Sample_Name': [], 'Sentrix_ID': [], 'Sentrix_Position': []}
# additional optional columns
addl_cols = []
if sample_type:
_dict['Sample_Type'] = []
addl_cols.append('Sample_Type')
if sample_sub_type:
_dict['Sample_Sub_Type'] = []
addl_cols.append('Sample_Sub_Type')
file_name_error_msg = "This .idat file does not have the right pattern to auto-generate a sample sheet: {0}"
for idat in idat_files:
# split string by '/', last element is local file name
try:
filename = str(idat).split("/")[-1]
split_filename = filename.split("_")
if split_filename[0].startswith('GSM'):
_dict['GSM_ID'].append(split_filename[0])
_dict['Sentrix_ID'].append(split_filename[1])
_dict['Sentrix_Position'].append(split_filename[2])
elif len(split_filename) == 3:
_dict['GSM_ID'].append("")
_dict['Sentrix_ID'].append(split_filename[0])
_dict['Sentrix_Position'].append(split_filename[1])
else:
raise ValueError(file_name_error_msg.format(idat))
except:
raise ValueError(file_name_error_msg.format(idat))
if sample_type:
_dict['Sample_Type'].append(sample_type)
if sample_sub_type:
_dict['Sample_Sub_Type'].append(sample_sub_type)
if matrix_file:
_dict['Sample_Name'] = sample_names_from_matrix(dir_path, _dict['GSM_ID'])
else:
# generate sample names
for i in range (1, len(_dict['GSM_ID']) + 1):
_dict['Sample_Name'].append("Sample_" + str(i))
df = pd.DataFrame(data=_dict)
df.to_csv(path_or_buf=(PurePath(dir_path, output_file)),index=False)
LOGGER.info(f"[!] Created sample sheet: {dir_path}/samplesheet.csv with {len(_dict['GSM_ID'])} GSM_IDs")
```
--------------------------------
### Run beta_bake via CLI
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/cli.html
Executes the beta_bake function to download and process a specific GEO dataset.
```bash
python -m methylprep beta_bake -i GSE74013 -d GSE74013
```
--------------------------------
### Initialize SigSet class
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/models/sigset.html
The constructor validates IDAT channel consistency and merges probe data with manifest dataframes.
```python
def __init__(self, sample, green_idat, red_idat, manifest, debug=False):
""" green_idat has .probe_means and .meta as main functions
and for extra info, use extra kwargs:
red= m.files.IdatDataset('9247377093_R02C01_Red.idat', m.models.Channel.RED, verbose=True, std_dev=True, nbeads=True)
"""
self.debug = debug
snps_read = {green_idat.n_snps_read, red_idat.n_snps_read}
if len(snps_read) > 1:
raise ValueError('IDAT files have a varying number of probes (comparing Grn to Red channel)')
if (str(green_idat.channel) != 'Grn' or str(red_idat.channel) != 'Red'):
raise ValueError("The IDAT files you supplied seem to be reversed. Check the order of your inputs to SigSet")
self.n_snps_read = snps_read.pop()
# DEBUG
#self.array_type = ArrayType.from_probe_count(self.n_snps_read)
# these next two should be unnecessary, because nothing should be reading idats downstream; use self.data_channel instead
#self.green_idat = green_idat
#self.red_idat = red_idat
self.data_channel = {'GREEN': green_idat.probe_means, 'RED': red_idat.probe_means} # indexed to illumina_ids
# illumina_ids are all II means, plus a stacked list of type-I-AddressA and type-I-AddressB means
self.sample = sample
self.man = manifest.data_frame # relevant columns are 'probe_type', AddressA_ID, AddressB_ID, index, Color_Channel
self.man = self.man[ ~self.man.index.str.startswith('rs') ] # snp_man covers these
self.snp_man = manifest.snp_data_frame.set_index('IlmnID')
self.ctl_man = manifest.control_data_frame
self.ctrl_green = self.ctl_man.merge(
green_idat.probe_means.astype('float32'),
how='inner', left_index=True, right_index=True)
self.ctrl_red = self.ctl_man.merge(
red_idat.probe_means.astype('float32'),
how='inner', left_index=True, right_index=True)
self.array_type = manifest.array_type
if self.array_type == ArrayType.ILLUMINA_MOUSE:
self.mouse_probes_mask = ( (self.man['design'] == 'Multi') | (self.man['design'] == 'Random') )
else:
self.mouse_probes_mask = None
self.address_code = {'AddressA_ID':'A', 'AddressB_ID':'B', 'A':'AddressA_ID', 'B':'AddressB_ID'}
"""
## SigSet EPIC
## - @IG probes: 49989 - 332 4145 70 7094 599 2958 ...
## - @IR probes: 92294 - 183 8040 1949 6152 833 89 ...
## - @II probes: 724612 - 6543 1596 3133 1011 3035 2837 ...
## - @oobG probes: 92294 - 138 277 107 218 232 80 ...
## - @oobR probes: 49989 - 1013 150 81 910 448 183 ...
## - @ctl probes: 635 ...
## - @pval: 866895 - 0.005141179 0.04914081 0.002757492 ...
SigSet 450k
@II 350076 ................... methylated 485512
@IG 46298 ... oobR 46298 ..... unmethylated 485512
@IR 89203 ... oobG 89203 ..... snp_methylated 65
.............................. snp_unmethylated 65
fg_green 396325 |vs| ibG 396374 (incl 40 + 9 SNPs) --(flattened)--> 442672
fg_red 439223 |vs| ibR 439279 (incl 40 + 16 SNPs) --(flattened)--> 528482
"""
if debug: print('DEBUG comparing [manifest probe_IDs vs idat probe_means]')
```
--------------------------------
### Partition Samples into Batches
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/_modules/methylprep/processing/pipeline.html
Organizes samples into batches based on a specified size, ensuring unique sample names and handling potential missing identifiers.
```python
batches = []
batch = []
sample_id_counter = 1
if batch_size:
if type(batch_size) != int or batch_size < 1:
raise ValueError('batch_size must be an integer greater than 0')
for sample in samples:
if sample_name and sample.name not in sample_name:
continue
# batch uses Sample_Name, so ensure these exist
if sample.name in (None,''):
sample.name = f'Sample_{sample_id_counter}'
sample_id_counter += 1
# and are unique.
if Counter((s.name for s in samples)).get(sample.name) > 1:
sample.name = f'{sample.name}_{sample_id_counter}'
sample_id_counter += 1
if len(batch) < batch_size:
batch.append(sample.name)
else:
batches.append(batch)
batch = []
batch.append(sample.name)
batches.append(batch)
else:
for sample in samples:
if sample_name and sample.name not in sample_name:
continue
# batch uses Sample_Name, so ensure these exist
if sample.name in (None,''):
sample.name = f'Sample_{sample_id_counter}'
sample_id_counter += 1
# and are unique.
if Counter((s.name for s in samples)).get(sample.name) > 1:
sample.name = f'{sample.name}_{sample_id_counter}'
sample_id_counter += 1
batch.append(sample.name)
batches.append(batch)
```
--------------------------------
### methylprep.run_series_list
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/source/modules.html
Downloads the IDATs and metadata for a list of series, creating metadata dictionaries and dataframes of sample beta_values.
```APIDOC
## methylprep.run_series_list
### Description
Downloads the IDATs and metadata for a list of series, creating metadata dictionaries and dataframes of sample beta_values.
### Parameters
- **list_file** (str) - Required - The name of the file containing a list of GEO_IDS and/or Array Express IDs.
- **path** (str) - Required - The path to the directory to download the data to.
- **dict_only** (bool) - Optional - Default False.
- **batch_size** (int) - Optional - Default 100.
```
--------------------------------
### Process IDAT files via command line
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/general_walkthrough.html
Executes the methylprep processing pipeline on a specified directory containing IDAT and manifest files.
```bash
>>> python -m methylprep process -d --all
```
--------------------------------
### Retrieve project sample sheet
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/IDE.html
Locates and parses the sample sheet file within a specified project directory. Returns a SampleSheet object containing the parsed information.
```python
from methylprep import get_sample_sheet
sample_sheet = get_sample_sheet(dir_path, filepath=None)
```
--------------------------------
### Run methylprep pipeline in a notebook
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/general_walkthrough.html
Executes the methylprep pipeline on a specified directory of raw data. Setting export=True is recommended to save output files for future use.
```python
from methylprep import run_pipeline
from pathlib import Path
filepath = Path('~/tutorial/GPL21145')
data_containers = run_pipeline(filepath, export=True)
```
--------------------------------
### methylprep.download.run_series_list
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/source/modules.html
Downloads IDATs and metadata for a list of series provided in a file.
```APIDOC
## methylprep.download.run_series_list(list_file, path, dict_only=False, batch_size=100, **kwargs)
### Description
Downloads the IDATs and metadata for a list of series, creating metadata dictionaries and dataframes of sample beta_values.
### Parameters
- **list_file** (str) - Required - The name of the file containing a list of GEO_IDS and/or Array Express IDs.
- **path** (str) - Required - The path to the directory to download the data to.
- **dict_only** (bool) - Optional - If True, only downloads data and creates dictionaries.
- **batch_size** (int) - Optional - Number of samples to process at a time (default 100).
```
--------------------------------
### Condensed GEO download and process workflow
Source: https://life-epigenetics-methylprep.readthedocs-hosted.com/en/latest/docs/faq.html
A complete workflow for downloading, filtering, and processing a GEO dataset.
```bash
python -m methylprep -v download -i GSE130030 -d GSE130030
# next, remove the treatment samples using `-c` and remove extra idats with `-s`
python -m methylprep -v meta_data -i GSE130030 -d GSE130030 --control -s
# finally, process it
python -m methylprep -v process -d GSE130030 --betas --m_value --no_export
```