NucProcess ---------- NucProcess is a Python program to perform single-cell Hi-C sequence processing of paired read data. This software takes paired FASTQ sequence read files (representing Hi-C data for only ONE cell), a reference genome sequence and knowledge of experimental parameters (restriction enzyme type and the range of size of DNA fragments sequenced in the library) to create processed, single-cell Hi-C contact files. By default, the output is generated in Nuc Chromatin Contact (NCC) format, which is a simple text based format described below. A processing report document and a contact map will also be generated for each run; in SVG format that can be viewed in most web browsers. To run NucProcess issue the 'nuc_process' command with the command line options described below. The options -i (input FASTQ files) and -g (genome reference) are mandatory, though its is usual to also use -re1 (primary restriction enzyme. Default is MboI), -o (root name of output files) and -re2 (secondary restriction enzyme in double-digest experiments). The 'nuc_contact_map' command takes the contact data from NCC format files to make all-chromosome contact map graphics in SVG format. This is automatically run on the main output of NucProcess, but can be run as required on any NCC format file. The 'nuc_contact_probability' command takes the contact data from one or more NCC format files to create log plots of contact probability versus sequence separation for intra chromosomal contacts. IMPORTANT: If the files for the genome index and the corresponding restriction enzyme (RE) cut location files have not been created they will be generated automatically by NucProcess so long as the -f option is specified. This option specifies the location of complete chromosomal sequences in FASTA format, typically by using a wild-card specification. Once these files are present the -f option need not be specified, but it will not trigger the re-creation of the files unless the -m or -x options are used (or if files are deleted). NOTE: the names of the chromosomes used by NucProcess are determined by the file names of the chromosome sequence files that were used to build the genome index and RE files (and these must match). The chromosome sequence files should be tagged with the chromosome name after an underscore but before the file extension. For example "chr7" will be taken from a file named "mm_ref_GRCm38.p2_chr7.fa". Citations --------- If you use NucProcess in published work, please cite the following reference: Stevens et al. Nature. 2017 Apr 6;544(7648):59-64 [PMID:28289288] Barcoded input -------------- The splitFastqBarcodes.py script is provided to split FASTQ files that represent many cells, each with a different barcode sequence, into separate paired read files. The script can be run as follows, specifying the names of the two paired, multiplex FASTQ read files after the script: python splitFastqBarcodes.py MULTIPLEXED_DATA_r_1.fq MULTIPLEXED_DATA_r_2.fq This will generate paired FASTQ files of the form: MULTIPLEXED_DATA_r_1_CGC.fq MULTIPLEXED_DATA_r_2_CGC.fq MULTIPLEXED_DATA_r_1_TAA.fq MULTIPLEXED_DATA_r_2_TAA.fq where the file names are tagged with the corresponding barcode sequence. These demultiplexed files can then be used as input to NucProcess, specifying only one barcode for each run. Running NucProcess ------------------ Typical first-time use: nuc_process -f /chromosomes/*.fa -o CELL_1 -v -a -k -re1 MboI -re2 AluI -s 150-2000 -n 12 -g /genome/GENOME_BUILD -i /data/SEQUENCING_DATA_r_?.fq Typical use thereafter: nuc_process -o CELL_1 -v -a -k -re1 MboI -re2 AluI -s 150-2000 -n 8 -g /genome/GENOME_BUILD -i /data/SEQUENCING_DATA_r_?.fq For the above commands: -f /chromosomes/*.fa states that all FASTA files (ending in .fa) in the /chromosomes/ directory will be used for creation of the genome index and RE cut site files -o specifies CELL_1 will be used for naming the output. In this case the main output contact file will be CELL_1.ncc -v specifies verbose output of processing progress -a specifies to generate ambiguous contact files: CELL_1_ambig.ncc in this case -k specifies to keep all the intermediate processing files: Filtered NCC files, clipped FASTQ files and the main Bowtie2 mapping SAM files -re1 is the primary restriction enzyme type at the ligation junction (see enzymes.conf) -re2 is the secondary restriction enzyme used to release the fragments (option not used for Nextera based protocol) -s is the valid molecule/fragment size range, as used in the DNA sequencing -n is the number of parallel CPU cores to use with Bowtie2 -g GENOME_BUILD is the root name for the Bowtie2 genome index without any file extension and in this case would refer to files GENOME_BUILD.1.bt2, GENOME_BUILD.rev.1.bt2 etc. -i SEQUENCING_DATA_r_?.fq is a wild-card expression matching the two input FASTQ files (though two separate file names, separated by a space can be specified). In this case the expression would match SEQUENCING_DATA_r_1.fq and SEQUENCING_DATA_r_2.fq - the paired sequence read files. To generate contact map graphics from an NCC format file: nuc_contact_map -i CELL_1.ncc This will generate the output graphics file CELL_1_contact_map.svg. However, the output file name maye be specified via the -o option. NCC data format --------------- The main NCC output format for contact data takes the form of space-separate lines, where each line represents a different chromosomal contact that pairs two chromosomal regions. This format specifies more than just the chromosome contact map: it also records the original read locations within the (ligated) primary RE digest fragments, strand information, ambiguity information and information to relate the data back to the original FASTQ input files. Accordingly this format can be used for data filtering and validation. The columns of NCC files correspond to: Name of chromosome A First base position of sequence read A Last base position of sequence read A 5' base position of primary RE fragment containing read A 3' base position of primary RE fragment containing read A The strand of sequence read A Name of chromosome B First base position of sequence read B Last base position of sequence read B 5' base position of primary RE fragment containing read B 3' base position of primary RE fragment containing read B The strand of sequence read B The number of the ambiguity group to which the paired reads belong The ID number of the read pair in the original FASTQ files Whether read pairs are swapped relative to original FASTQ files i.e: chr_a start_a end_a re1_a_start re1_a_end strand_a chr_b start_b end_b re1_b_start re1_b_end strand_b ambig_group pair_id swap_pair For example two lines could be: chr10 100002828 100002899 100002733 100003107 + chr10 100015771 100015700 100015676 100016001 - 1498612 1534464 0 chr10 100007729 100007658 100007551 100008354 - chr13 107185630 107185700 107184984 107185698 + 1009602 1032357 1 Command line options for nuc_process ------------------------------------ usage: nuc_process [-h] [-i FASTQ_FILE [FASTQ_FILE ...]] [-g GENOME_FILE] [-g2 GENOME_FILE_2] [-re1 ENZYME] [-re2 ENZYME] [-s SIZE_RANGE] [-n CPU_COUNT] [-r COUNT] [-o NCC_FILE] [-oa NCC_FILE] [-or REPORT_FILE] [-b EXE_FILE] [-q SCHEME] [-qm MIN_QUALITY] [-m] [-p] [-pt PAIRED_READ_TAGS PAIRED_READ_TAGS] [-x] [-f FASTA_FILES [FASTA_FILES ...]] [-f2 FASTA_FILES_2 [FASTA_FILES_2 ...]] [-a] [-k] [-sam] [-l SEQUENCE] [-z] [-v] [-hc HOM_CHROMO_TSV_FILE] [-u] [-c GENOME_COPIES] Chromatin contact paired-read Hi-C processing module for Nuc3D and NucTools optional arguments: -h, --help show this help message and exit -i FASTQ_FILE [FASTQ_FILE ...] Input paired-read FASTQ files to process. Accepts wildcards that match paired files. If more than two files are input, processing will be run in batch mode using the same parameters. -g GENOME_FILE Location of genome index files to map sequence reads to without any file extensions like ".1.b2" etc. A new index will be created with the name if the index is missing and genome FASTA files are specified -g2 GENOME_FILE_2 Location of secondary genome index files for hybrid genomes. A new index will be created with the name if the index is missing and genome FASTA files are specified -re1 ENZYME Primary restriction enzyme (for ligation junctions). Default: MboI. Available: AluI, BglII, DpnII, HindIII, MboI -re2 ENZYME Secondary restriction enzyme (if used). Available: AluI, BglII, DpnII, HindIII, MboI -s SIZE_RANGE Allowed range of sequenced molecule sizes, e.g. "150-1000", "100,800" or "200" (no maximum) -n CPU_COUNT Number of CPU cores to use in parallel -r COUNT Minimum number of sequencing repeats required to support a contact -o NCC_FILE Optional output name for NCC format chromosome contact file. This option will be ignored if more than two paired FASTA files are input (i.e. for batch mode); automated naming will be used instead. -oa NCC_FILE Optional output name for ambiguous contact NCC file. This option will be ignored if more than two paired FASTA files are input (i.e. for batch mode); automated naming will be used instead. -or REPORT_FILE Optional output name for SVG format report file. This option will be ignored if more than two paired FASTA files are input (i.e. for batch mode); automated naming will be used instead. -b EXE_FILE Path to bowtie2 (read aligner) executable (will be searched for if not specified) -q SCHEME Use a specific FASTQ quality scheme (normally not set and deduced automatically). Available: phred33, phred64, solexa -qm MIN_QUALITY Minimum acceptable FASTQ quality score in range 0-40 for clipping 3' end of reads. Default: 10 -m Force a re-mapping of genome restriction enzyme sites (otherwise cached values will be used if present) -p The input data is multi-cell/population Hi-C; single- cell processing steps are avoided -pt PAIRED_READ_TAGS PAIRED_READ_TAGS When more than two FASTQ files are input (batch mode), the subtrings/tags which differ between paired FASTQ file paths. Default: r_1 r_2 -x, --reindex Force a re-indexing of the genome (given appropriate FASTA files) -f FASTA_FILES [FASTA_FILES ...] Specify genome FASTA files for genome index building (accepts wildcards) -f2 FASTA_FILES_2 [FASTA_FILES_2 ...] A second set of genome FASTA files for building a second genome index when using hybrid strain cells (accepts wildcards). -a Whether to report ambiguously mapped contacts -k Keep any intermediate files (e.g. clipped FASTQ etc). -sam Write paired contacts files to SAM format -l SEQUENCE Seek a specific ligation junction sequence (otherwise this is guessed from the primary restriction enzyme) -z GZIP compress any output FASTQ files -v, --verbose Display verbose messages to report progress -hc HOM_CHROMO_TSV_FILE, --homologous_chromos HOM_CHROMO_TSV_FILE File path specifying whitespace-separated pairs of homologous chromosome names (to match first word in header lines of genome sequence FASTQ files) and corresponding pair name (e.g. "chrX"). Required for hybrid strain analysis. See -g2 option. -u Whether to only accept uniquely mapping genome positions and not attempt to resolve certain classes of ambiguous mapping where a single perfect match is found. -c GENOME_COPIES Number of whole-genome copies, e.g. for S2 phase; Default 1 unless homologous chromosomes (-hc) are specified for hybrid genomes. Note enzymes.conf can be edited to add further restriction enzyme cut-site definitions. Command line options for nuc_contact_map ---------------------------------------- usage: nuc_contact_map [-h] [-i NCC_FILE [NCC_FILE ...]] [-o SVG_FILE_TAG] [-w SVG_WIDTH] [-s BIN_SIZE] [-b] [-c RGB_COLOR] Chromatin contact (NCC format) Hi-C contact map display module for Nuc3D and NucTools optional arguments: -h, --help show this help message and exit -i NCC_FILES Input NCC format chromatin contact file(s). Wildcards accepted -o SVG_FILE Optional name tag to put at end of SVG format contact map file. Use "-" to print SVG to stdout rather than make a file. Default: "_contact_map" -w SVG_WIDTH SVG document width -s BIN_SIZE Sequence region size represented by each small square (the resolution) in megabases. Default is 5 kb -b Specifies that the contact map should have a black background (default is white) -c RGB_COLOR Optional main color for the contact points as a 24-bit hexidecimal RBG code e.g. "#0080FF" (with quotes) For further help email tjs23@cam.ac.uk or wb104@cam.ac.uk Command line options for nuc_contact_probability ------------------------------------------------ usage: nuc_contact_probability [-h] [-i NCC_FILE [NCC_FILE ...]] [-o SVG_FILE] [-w SVG_WIDTH] [-s KB_BIN_SIZE] Chromatin contact (NCC format) probability vs sequence separation graph module for Nuc3D and NucTools optional arguments: -h, --help show this help message and exit -i NCC_FILE [NCC_FILE ...] Input NCC format chromatin contact file(s). Wildcards accepted -o SVG_FILE Output SVF format file. Use "-" to print SVG to stdout rather than make a file. -w SVG_WIDTH SVG document width -s KB_BIN_SIZE Sequence region size in kilobases for calculation of contact probabilities. Default is 100 (kb)
kwohlfahrt/nuc_processing
Chromatin contact paired-read single-cell Hi-C processing module for Nuc3D and NucTools
PythonLGPL-3.0