-
Notifications
You must be signed in to change notification settings - Fork 2
Home
To combine barcodes with 10X we insert transcribable lineage barcodes into cells and then perform 10X on those barcoded cells. 10X will capture the transcribed lineage barcode and attach it's cell barcode to it. From the completed full length 10X library prep we use primers to amplify the barcode sequences out of the 10X material keeping both the 10X cell barcode and our lineage barcode in the amplified product. We then use paired end sequencing to capture the 10X barcode (Read1) as well as our barcode (Read2). In the end we have sequencing runs containing the gene expression data of the cells and a sequencing run with all the barcodes. This wiki walks you through how to use cellranger and this code to combine these sequencing runs so that you can attribute lineage barcodes to specific cells sequenced in your 10X run. Specifically the scripts provided here take the barcode sequencing run and modifies the R2 files as well as creates a file with all detected barcodes so that we can run the Feature Barcode Analysis (https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/3.1/using/feature-bc-analysis) pipeline written by 10X. Additionally if a portion of barcoded cells used in the 10x were used in other experiments where gDNA was collected this pipeline allows you do identify the same barcodes found in the 10x libraries in gDNA libraries.
All this code is run on the pmacs cluster as starcode is very memory intensive. NOTE: If you have already installed this file on the cluster and compiled starcode start at step 9
- If you don’t already have it I recommend you download cyberduck (https://cyberduck.io) to easily interface with the cluster
- Once cyberduck is installed, open it, hit open connection, in the top drop-down menu select "SFTP(SSH File Transfer Protocol)", for server type in "mercury.pmacs.upenn.edu", enter your PMACS username and password, and hit connect.
- Download all the files in this GitHub repository by hitting the green "Code" button above on this page, and select "Download ZIP"
- Drag and drop the downloaded file "BarcodeAnalysis-main" into your home directory in cyberduck (for future steps to work without modification this needs to be in your home directory, in cyberduck the drop down should read "/home/username")
- Change the file name from "BarcodeAnalysis-main" to simply "BarcodeAnalysis"
- To set up the virtual environment containing the necessary dependencies, open up terminal and type in
ssh <username>@consign.pmacs.upenn.edureplacing<username>with your PMACS username, and hit enter. Then enter your PMACS password and hit enter - You should now be in a session on the terminal.Enter
bsub -Is bashto start up a node - The first step is to compile starcode by running:
make -C /home/<username>/BarcodeAnalysis/starcodereplacing<username>with your PMACS username - We will now make sure we are in the correct python version by entering the command
module load python/3.6.3. If this worked you should be able to enterpython --versionand see the outputPython 3.6.3 - Enter
cd /home/<username>/BarcodeAnalysis/replacing<username>with your PMACS username. (tip: to get paths of folders I usually right click on them in cyberduck, select "copy URL", and select either choice, and keep everything after ".edu") - Now set up your virtual environment by entering the line
python -m virtualenv bcEnv - Next activate your virtual environment by entering
source /home/<username>/BarcodeAnalysis/bcEnv/bin/activate, replacing with your PMACS username - We can now install the necessary packages in this environment by running:
pip install -r /home/<username>/BarcodeAnalysis/ExtractBarcodes/requirements.txt, replacing<username>with your PMACS username - Your virtual environment is now ready. If you ever want to leave the virtual environment simply enter
deactivate, every time you want to active the environment enter:source /home/<username>/BarcodeAnalysis/bcEnv/bin/activate - You also need to be able to run CellRanger on the cluster, you can access cellranger by either installing it as described here:https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/latest/installation , or you can use the version of cellranger on the cluster. To use the version on the cluster, simply log in, activate a node by typing in
bsub -Is bash, hitting enter, and then enteringmodule load cellranger/5.0.1(you can see if there are other versions by typingmodule avail). If you are downloading cell ranger make sure you add it you your path every time before running a cellranger command by entering this commandexport PATH=/path/to/cellranger-3.1.0:$PATHreplacing with your path to the cellranger folder. If you are using the module on the cluster then just runmodule load cellranger/5.0.1every time before running cellranger. Note: I downloaded cell ranger so my example files will useexport PATH=/path/to/cellranger-3.1.0:$PATH.
Here I walk through how to generate the fastq files from the raw sequencing data collected by sequencing the 10X library, and your lineage barcode library. If you are confident you have properly generated your fastq files then you can skip to the "Modify R2 in the Lineage Barcode FASTQ Files" step.
- The first step is to download your data from BaseSpace on to your computer or a hard drive for safe keeping. You can do this by running the following command
bs download run -i <run ID> -o /path/to/download/location(you need to have the base space command line interface installed to do this)you can find the run ID in the URL of the run on base space. For example the run ID in this URL: "https://basespace.illumina.com/run/190309133/10Xbarcode_20190808/details" is 190309133. You will have to run this command for all the sequencing runs your performed on your 10x library and lineage barcode library. - Once you have downloaded the sequencing runs you can copy the data to the cluster where we will process it (Note: it is possible to download you data straight to the cluster but I like having a local copy of the raw data). To do this run the command
rsync -rav /Path/to/downloaded/seq/data <username>@mercury.pmacs.upenn.edu:/Path/to/location/on/clusterin terminal replacing the paths to your specific path and replacing with your PMACS username. Hit enter and enter your PMACS password and hit enter again. - Once the data is on the cluster we can use CellRanger to make fastq files. The instructions on how to do so can be found here: https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/latest/using/mkfastq. I will also detail how I made my fastqs in this document.
- To run cellranger mkfastq you need two files for every sequencing run you did. The first file is a mkfastq.sh file which contain the command to run cellranger mkfastq and the second is a sample sheet (samplesheet.csv) which contains the name of the samples and their respective barcodes. I have an example of both these files in this repository which you can find in the "ExtractBarcodes/ExampleFiles" folder. More details about these files and how to make them can be found here: https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/latest/using/mkfastq (NOTE: my files were were made to work with version 3.1 of cellranger so you will need to make modification for the newer versions of cellranger). The format for all the 10x library sequencing files will all be very similar, but it will be different for the barcode sequencing library since they do not use the 10x indexes. One of the key differences is that you actually have to input the index sequences, which should both (i5 and i7) be reverse complemented in this file. Examples for the files to use when running mkfastq on the lineage barcode library can also be found in the "ExtractBarcodes/ExampleFiles" folder and are called BC_mkfastq.sh and BC_samplesheet.csv.
- Once you have made all the files you can run them by first activating a node on the cluster by running
bsub -Is bashand then run this line for each sequencing runbsub < /path/to/mkfastq.sh. These should run for a while as they generate fastqs.
-
For the gDNA we can't use cellranger to make our fastqs and instead need to use bcl2fastq. To do this we need the bcl2fastq.sh file which contains the command to run bcl2fastq and gDNA_BC_samplesheet.csv which contains the names of the samples and their respective indices (you can find examples of this file in the directory: "ExtractBarcodes/ExampleFiles". In the bcl2fastq you will need to change paths for where the inputs are and where you want the fastqs to be output. In the sample sheet enter the read lengths (usually 151 for gDNA barcodes) in the [Reads] section. In the [Data] section under Sample_ID enter all the sample names (I recommend naming the gDNA_) as well as the reverse complement of the i7 in the column named "index" and the revers complement of the i5 in the column named index2.
-
To run bcl2fast enter
bsub < /path/to/bcl2fastq.sh
NOTE: you only need to run the steps outlined in this section once per sample. If you already have these results skip to the "Run Barcode Extraction" section below.
- To make sure we extract the most confidently identified barcodes we start by generating some plots to help guide what inputs to use when extracting the barcodes. All paths and variables than can be changed are located in the paths_and_variables.json file (located here: /home//BarcodeAnalysis/ExtractBarcodes/paths_and_variables.json) 2.The first step to generating these metrics is to fill out this paths_and_variables.json file. The inputs are mostly explained in this file and give an example path of where things are likely located. Here is the inputs asked for in the order they appear in the file:
starcode_path: The path to the starcode folder
Fastqfolder10x: The path to the folder containing all the different 10x sample folders and their FASTQs (if FASTQs were generated will cellranger the path should look something like:"/project/shafferslab/Guillaume/10X_exp1_reanalysis/20190808_10X1_BC_r1_r2_seq1_fastq/outs/fastq_path/10x")
FastqfoldergDNA: this is the same as above but giving the path to the gDNA fastqs instead of the 10X fastqs. If you don't have gDNA data input ""
Outfolder: folder where all outputs for extracting barcodes are saved
barcodeSource: Input if you are analyzing only 10x barcodes ("10x"), only gDNA barcodes ("gDNA"), or both "both"
strtseq: This is the sequence that is common before all barcodes. you should not have to change this if you are using the original barcodes designed by Ben Emert.
strtseq_revcomp: This is just the reverse complement of strtseq entered above
bclen: Ignore for now, leave input as 70
sc_mm: Ignore for now, leave input as 8
startseqMatch: ignore for now, leave input as 70
num_of_barcodes_to_use: this is the number of barcodes used to calculate Levantine distances depending on the metrics you use. we suggest using an input of ~500 doing may slightly increase the accuracy of the test results, but will increase the length of time it will take to generate the metric plots.
GSAMP: Define which samples should be run together in starcode (any samples that could contain cells from the same lineage should be grouped together in a list within the GSAMP list(see example in the paths_and_variables.json file). NOTE: names used here should match the names that will be assinged to samples in line 229 and 250 in the step3.py script (should be everything before lane (_L0)).
spike_in_added: Does your sample contain spike-ins "yes" or "no". If "no" ignore the rest of the variables:'spike_in_seqs','spike_in_value' and 'spike_in_color'. If "yes" fill out 'spike_in_seqs','spike_in_value' and 'spike_in_color'. The order of 'spike_in_seqs' has to match 'spike_in_value' and what ever color you picked for 'spike_in_color'.
spike_in_seqs: A list of known spike-in sequences without the startseq
spike_in_value: A list with the number of cells added for each sequence in 'spike_in_seqs'.
spike_in_color: A list with a specific color for each sequence in 'spike_in_seqs'.
-
Make sure your are in the environment bcEnv (if not activate it by running:
source /home/<username>/BarcodeAnalysis/bcEnv/bin/activatewith your username) and on a sub node (run:bsub -Is bashif not) and that you are currently in the directory containing the paths_and_variables.json file you want to use. For example, if you did not move the paths_and_variables.json file you need to be in the directorycd /home/<username>/BarcodeAnalysis/ExtractBarcodes/. To run the first step enter the linebsub -e step1.e -o step1.o -M 256000 python /home/<username>/BarcodeAnalysis/ExtractBarcodes/Scripts/Step1_optional.pymaking sure to substitute in your user name for . -
Once the above is done running run the next step by entering the line
bsub -e step2.e -o step2.o -M 256000 python /home/<username>/BarcodeAnalysis/ExtractBarcodes/Scripts/Step2_optional.pymaking sure to substitute in your user name for
-
Using the plots generated above we can now decide on the final values we want to use to extract barcodes. NOTE: you only need to run the steps outlined in "Get barcode metrics" once per sample.
-
The first value to decide on is to select how much of the barcode you want to use in your analysis. This is the "bclen" parameter in the paths_and_variables.json file. To do this look at the LV_distance plots generated in the Get barcode metrics section (these can be found in the following path: /output/path/specified/in/paths_and_variables.json/check_barcodes/LV_distance). Look for what length allows you to easily see a difference between barcodes that are the same and those that are different. This means that you should be able to see a small distribution that has very small levenshtein distances and then a bigger one with large levenshtein distances. Usually this will be the largest barcode size (70) unless there are significant read quality issues at the end of the sequence. Input what ever you decide for bclen in the paths_and_variables.json file.
-
The next parameter to decide is how many mistakes to allow between barcodes and still call them the same barcode. This is the "sc_mm" parameter in the paths_and_variables.json file. To figure out the value of this parameter look at the same LC_distance histogram we used in step 2 above, but only look at the ones that correspond the the barcode length you selected. The value you want to select here is the highest Levenshtein distance in the first small distribution (the highest value accepted by starcode is 8). Input what ever you decide for sc_mm in the paths_and_variables.json file.
-
The last parameter to decide on is how many mismatches in the start sequence (the sequence that identifies the start of a barcode). This is the "startseqMatch" parameter in the paths_and_variables.json file. To figure this out look at the Match_to_StartSeq plots generated in the Get barcode metrics section (these can be found in the following path: /output/path/specified/in/paths_and_variables.json/check_barcodes/Match_to_StartSeq). Based on this plot select the percent match that you believe is accurately detecting start sequences without sacrificing too much data (there is no clear way to select the best value). Input what ever you decide for startseqMatch in the paths_and_variables.json file.
-
Make sure you have properly saved your selected values in the paths_and_variables.json file before moving to the next step.
-
We can now extract the barcode by running script3.py. Make sure your are in the environment bcEnv (if not activate it by running:
source /home/<username>/BarcodeAnalysis/bcEnv/bin/activatewith your username) and on a sub node (run:bsub -Is bashif not) and that you are currently in the directory containing the paths_and_variables.json file you want to use. For example, if you did not move the paths_and_variables.json file you need to be in the directorycd /home/<username>/BarcodeAnalysis/ExtractBarcodes/. To run the first step enter the linebsub -e step3.e -o step3.o -M 256000 python /home/<username>/BarcodeAnalysis/ExtractBarcodes/Scripts/Step3.pymaking sure to substitute in your user name for . -
-Optional- Once Step3.py is done running, we can analyze the barcodes collapsed by starcode between samples in GSAMP. Step4 will create a FeatureReference_filtered_group.csv file for every group in GSAMP and will also plot comparisons between the samples. Run Step4 by entering the line
python /home/<username>/BarcodeAnalysis/ExtractBarcodes/Scripts/Step4_optional.pymaking sure to substitute in your user name for . -
Once done you should find the outputs needed to run cell ranger in the folder: /output/path/specified/in/paths_and_variables.json/check_barcodes/CellRanger_inputs. We suggest using the FeatureReference_filtered.csv when running starcode as this removes barcodes with a lot of base repeats, which shouldn't exist in the barcode library, and there for are unlikely to be real barcodes.
To generate count matrices form this data we now need to run cellranger count with feature barcodes as described here: https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/latest/using/feature-bc-analysis.
- The first step is to create a library.csv for each sample. The first columns of this file contains the paths to the fastq files of each run that contained a given samples, the second column has the name of the sample, the last column indicates if that row are reads for "Gene Expression" or if they are from the lineage barcodes in which case you enter "Custom". An example of this file is in this repository in the "Extract10xBarcodes/ExampleFiles" folder.
- Next create a count.sh file for each sample that contains the command to run the count. An example of this file is in this repository in the "Extract10xBarcodes/ExampleFiles" folder.
- Run count for each sample by running
bsub < /path/to/count.shon an active node on the cluster for each count.sh file. This should result in having count matrices for each sample.
If you have different samples that you would like to put into one count matrix, such as replicates, then you can aggregate your count matrices using cellranger aggr as described here: https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/6.0/what-is-cell-ranger
- To run cellranger aggr you must first make an aggregate csv that contains the sample name in the first column, the path to the molecule_info.h5 file (generated in the count step above) for that sample in the second column, and any additional data you might want about the sample in additional columns (anything after the second column is optional). You will also need a file to run the aggregate function on the cluster called aggregate.sh. An example of both these files is in this repository in the "Extract10xBarcodes/ExampleFiles" folder.
- Next just run cellranger aggr using the command
bsub < /path/to/aggregate.shon an active node in the cluster. You should now have a single count matrix containing the cells from all your samples ready for further custom analysis in R or python.