HackMD
    • Sharing Link copied
    • /edit
    • View mode
      • Edit mode
      • View mode
      • Book mode
      • Slide mode
      Edit mode View mode Book mode Slide mode
    • Customize slides
    • Note Permission
    • Read
      • Owners
      • Signed-in users
      • Everyone
      Owners Signed-in users Everyone
    • Write
      • Owners
      • Signed-in users
      • Everyone
      Owners Signed-in users Everyone
    • Commenting & Invitee
    • Publishing
      Please check the box to agree to the Community Guidelines.
      Everyone on the web can find and read all notes of this public team.
      After the note is published, everyone on the web can find and read this note.
      See all published notes on profile page.
    • Commenting Enable
      Disabled Forbidden Owners Signed-in users Everyone
    • Permission
      • Forbidden
      • Owners
      • Signed-in users
      • Everyone
    • Invitee
    • No invitee
    • Options
    • Versions and GitHub Sync
    • Transfer ownership
    • Delete this note
    • Note settings
    • Template
    • Insert from template
    • Export
    • Dropbox
    • Google Drive Export to Google Drive
    • Gist
    • Import
    • Dropbox
    • Google Drive Import from Google Drive
    • Gist
    • Clipboard
    • Download
    • Markdown
    • HTML
    • Raw HTML
Menu Note settings Sharing Help
Menu
Options
Versions and GitHub Sync Transfer ownership Delete this note
Export
Dropbox Google Drive Export to Google Drive Gist
Import
Dropbox Google Drive Import from Google Drive Gist Clipboard
Download
Markdown HTML Raw HTML
Back
Sharing
Sharing Link copied
/edit
View mode
  • Edit mode
  • View mode
  • Book mode
  • Slide mode
Edit mode View mode Book mode Slide mode
Customize slides
Note Permission
Read
Owners
  • Owners
  • Signed-in users
  • Everyone
Owners Signed-in users Everyone
Write
Owners
  • Owners
  • Signed-in users
  • Everyone
Owners Signed-in users Everyone
Comment & Invitee
Publishing
Please check the box to agree to the Community Guidelines.
Everyone on the web can find and read all notes of this public team.
After the note is published, everyone on the web can find and read this note.
See all published notes on profile page.
More (Comment, Invitee)
Commenting Enable
Disabled Forbidden Owners Signed-in users Everyone
Permission
Owners
  • Forbidden
  • Owners
  • Signed-in users
  • Everyone
Invitee
No invitee
   owned this note    owned this note      
Published Linked with GitHub
Like BookmarkBookmarked
Subscribed
  • Any changes
    Be notified of any changes
  • Mention me
    Be notified of mention me
  • Unsubscribe
Subscribe
# Map/align metagenome short reads to a reference database BBTools is a powerful collection of commands to analyze metagenomic sequences. Download the software [here](https://sourceforge.net/projects/bbmap/) and find a user guide [here](https://jgi.doe.gov/data-and-tools/bbtools/bb-tools-user-guide/). We will use the tool BBMap to align short reads to a reference database. The idea is to retrieve reads that belong to certain microbial clade thereby removing reads of organisms that are not of interest. In this example we are interested in sulfate reducers. We have already created a [reference database](https://hackmd.io/4V9elLP3Rr2gUAHLf2DPuw) with genomes belonging to known sulfate-reducing bacteria. ## Map reads to reference database To map metagenomic reads to a reference database we first build an index. You created your 'reads.fna' in the chapter [Create genome reference database](https://hackmd.io/9JxrczlwRtCm_hFPtyh1hw?view). ``` module load bbtools module load samtools bbmap.sh ref=your-database-file.fna ``` Note:`ref=your-database-file.fna` is named differently and might have the extension 'gz', e.g. `ref=desulfo_genomes_DB.fna.gz`. The output looks like this: ``` No output file. Writing reference. Executing dna.FastaToChromArrays2 [desulfo_genomes_DB.fna, 1, writeinthread=false, genscaffoldinfo=true, retain, waitforwriting=false, gz=true, maxlen=536670912, writechroms=true, minscaf=1, midpad=300, startpad=8000, stoppad=8000, nodisk=false] Set genScaffoldInfo=true Writing chunk 1 Writing chunk 2 Set genome to 1 Loaded Reference: 0.154 seconds. Loading index for chunk 1-2, build 1 No index available; generating from reference genome: /users/eruff/Spart2020_Student/genome_databases/desulfo/ref/index/1/chr1-2_index_k13_c2_b1.block Indexing threads started for block 0-2 Indexing threads finished for block 0-2 Generated Index: 177.753 seconds. Finished Writing: 30.043 seconds. No reads to process; quitting. Total time: 224.825 seconds. ``` The command produces a directory called `ref` with several sub-directories. Now we can map the metagenomic reads from metagenomes to the index using the following command. This will create two `bam` files (short for Binary Alignment Map), one containing mapped and one containing unmapped reads. Note: The mapping of 1GB metagenomic reads takes around 10 minutes on an average cluster using 12 threads. ``` bbmap.sh in=/groups/spartina2020/genome_databases/MG.fastq outu=MG_unmapped.bam outm=MG_mapped.bam ``` In case the metagenome consists of paired reads it is recommended to process them as such. Most bbtools can handle paired reads, check [here](https://jgi.doe.gov/data-and-tools/bbtools/bb-tools-user-guide/usage-guide/). Paired reads can be mapped and turned into a single interleaved file using this command: ``` bbmap.sh in1=/groups/spartina2020/genome_databases/MG_R1.fastq in2=/groups/spartina2020/genome_databases/MG_R2.fastq outu=MG_unmapped.bam outm=MG_mapped.bam ``` The console output looks somewhat like this: ``` ------------------ Results ------------------ Genome: 1 Key Length: 13 Max Indel: 16000 Minimum Score Ratio: 0.56 Mapping Mode: normal Reads Used: 7261066 (1063955727 bases) Mapping: 602.129 seconds. Reads/sec: 12058.99 kBases/sec: 1766.99 ... Total time: 639.210 seconds. ``` The bam file of the mapped reads is converted back to fastq files using samtools ``` samtools fastq MG_mapped.bam > MG_mapped.fastq ``` ## Process multiple sequences at once Because it is inefficient to process every sample on its own, let's tell bbmap to map/align multiple metagenomes on one index file. For this we will create a `for loop`. Note: Even if multiple input files are processed the index file needs to be created only once. Find more info [here](https://github.com/BioInfoTools/BBMap/blob/master/sh/bbmap.sh). ### Create a list with all metagenomes that need to be processed This can be done using a number of commands. In ths case we just need a list of the names of all `fastq.gz` files. This will do: `ls /path/to/dir/*.fastq.gz > /path/to/dir/file_list.txt` If you need a list with files and their path names use, e.g this: ``` find /path/to/dir/ -type f -name '*.fastq.gz' -exec cat {} + > /path/to/dir/file_list.txt ``` ### Now let's write a loop The following `for loop` will take each metagenome in the above list, and maps/aligns it to our index. (Remeber, we do not have to specify the location of the index as long as we are in the same directory as the index). The loop creates separate bam files for the `mapped` and `unmapped` reads. The bam files will be named after the metagenome they are coming from so that we don't loose track. Find more info on a very similar loop [here](https://www.biostars.org/p/381398/). ``` for file in $(cat IGERT_mgm_list.txt); do bbmap.sh in=/groups/spart2020/maptest/IGERT_mgm/${file} outm=./bam/${file}_mapped.bam outu=./bam/${file}_unmapped.bam; done ```

Import from clipboard

Advanced permission required

Your current role can only read. Ask the system administrator to acquire write and comment permission.

This team is disabled

Sorry, this team is disabled. You can't edit this note.

This note is locked

Sorry, only owner can edit this note.

Reach the limit

Sorry, you've reached the max length this note can be.
Please reduce the content or divide it to more notes, thank you!

Import from Gist

Import from Snippet

or

Export to Snippet

Are you sure?

Do you really want to delete this note?
All users will lost their connection.

Create a note from template

Create a note from template

Oops...
This template is not available.


Upgrade

All
  • All
  • Team
No template found.

Create custom template


Upgrade

Delete template

Do you really want to delete this template?

This page need refresh

You have an incompatible client version.
Refresh to update.
New version available!
See releases notes here
Refresh to enjoy new features.
Your user state has changed.
Refresh to load new user state.

Sign in

Forgot password

or

By clicking below, you agree to our terms of service.

Sign in via Facebook Sign in via Twitter Sign in via GitHub Sign in via Dropbox

New to HackMD? Sign up

Help

  • English
  • 中文
  • Français
  • Deutsch
  • 日本語
  • Español
  • Català
  • Ελληνικά
  • Português
  • italiano
  • Türkçe
  • Русский
  • Nederlands
  • hrvatski jezik
  • język polski
  • Українська
  • हिन्दी
  • svenska
  • Esperanto
  • dansk

Documents

Tutorials

Book Mode Tutorial

Slide Mode Tutorial

YAML Metadata

Contacts

Facebook

Twitter

Discord

Feedback

Send us email

Resources

Releases

Pricing

Blog

Policy

Terms

Privacy

Cheatsheet

Syntax Example Reference
# Header Header 基本排版
- Unordered List
  • Unordered List
1. Ordered List
  1. Ordered List
- [ ] Todo List
  • Todo List
> Blockquote
Blockquote
**Bold font** Bold font
*Italics font* Italics font
~~Strikethrough~~ Strikethrough
19^th^ 19th
H~2~O H2O
++Inserted text++ Inserted text
==Marked text== Marked text
[link text](https:// "title") Link
![image alt](https:// "title") Image
`Code` Code 在筆記中貼入程式碼
```javascript
var i = 0;
```
var i = 0;
:smile: :smile: Emoji list
{%youtube youtube_id %} Externals
$L^aT_eX$ LaTeX
:::info
This is a alert area.
:::

This is a alert area.

Versions

Versions and GitHub Sync

Sign in to link this note to GitHub Learn more
This note is not linked with GitHub Learn more
 
Add badge Pull Push GitHub Link Settings
Upgrade now

Version named by    

More Less
  • Edit
  • Delete

Note content is identical to the latest version.
Compare with
    Choose a version
    No search result
    Version not found

Feedback

Submission failed, please try again

Thanks for your support.

On a scale of 0-10, how likely is it that you would recommend HackMD to your friends, family or business associates?

Please give us some advice and help us improve HackMD.

 

Thanks for your feedback

Remove version name

Do you want to remove this version name and description?

Transfer ownership

Transfer to
    Warning: is a public team. If you transfer note to this team, everyone on the web can find and read this note.

      Link with GitHub

      Please authorize HackMD on GitHub

      Please sign in to GitHub and install the HackMD app on your GitHub repo. Learn more

       Sign in to GitHub

      HackMD links with GitHub through a GitHub App. You can choose which repo to install our App.

      Push the note to GitHub Push to GitHub Pull a file from GitHub

        Authorize again
       

      Choose which file to push to

      Select repo
      Refresh Authorize more repos
      Select branch
      Select file
      Select branch
      Choose version(s) to push
      • Save a new version and push
      • Choose from existing versions
      Available push count

      Upgrade

      Pull from GitHub

       
      File from GitHub
      File from HackMD

      GitHub Link Settings

      File linked

      Linked by
      File path
      Last synced branch
      Available push count

      Upgrade

      Danger Zone

      Unlink
      You will no longer receive notification when GitHub file changes after unlink.

      Syncing

      Push failed

      Push successfully