run_cutadapt: Run Cutadapt

Description Usage Arguments Value Examples

View source: R/run_cutadapt.R

Description

Run the Cutadapt tool to remove sequencing adapters and low quality bases.

Usage

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
run_cutadapt(
  mate1 = NULL,
  mate2 = NULL,
  mate1.out = NULL,
  mate2.out = NULL,
  quality = NULL,
  nextseq = FALSE,
  minimum = NULL,
  trim.only = FALSE,
  cut = NULL,
  adapter1 = NULL,
  adapter2 = NULL,
  polyA = NULL,
  parallel = FALSE,
  cores = 4,
  cutadapt = NULL,
  version = FALSE
)

Arguments

mate1

List of the paths to files containing to the forward reads

mate2

List of the paths to files containing to the reverse reads

mate1.out

List of paths to the files to write the trimmed forward reads

mate2.out

List of paths to the files to write the trimmed reverse reads

quality

The lower limit for the phred score

nextseq

Was the sequence data generated on a NextSeq 500, trims dark cycle bases appearing as high-quality G bases

minimum

The length at which a trimmed read will be discarded

trim.only

Only keep reads that have had adapters trimmed

cut

Remove the first 'n' bases form the 5' end of the forward read

adapter1

Sequence for the adapter for the forward read

adapter2

Sequence for the adapter for the reverse read

polyA

Number of A's

parallel

Run in parallel, default set to FALSE

cores

Number of cores/threads to use for parallel processing, default set to 4

cutadapt

Path to the Cutadapt program, required

version

Returns the version number

Value

A file with the Cutadapt commands and creates a directory of adapter and quality trimmed reads

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
 ## Not run: 
 # Version number
 run_cutadapt(cutadapt = cutadapt.path,
                         version = TRUE)

# Trimmed reads directory
trimmed.reads.dir <- "trimmed_reads"
#Create the directory for the trimmed reads
dir.create(trimmed.reads.dir, showWarnings = FALSE)

read1.pattern <- "*_R1_001.fastq.gz$"
read2.pattern <- "*_R2_001.fastq.gz$"

reads.path <- "/export/buzz2/gpfrawdata/nextseq01/FastQ/2019/NeilBulleid/
               1466_NS205_0325_MarieAnnPringle_20190313"

mate1 <- list.files(path = reads.path,
                    pattern = read1.pattern,
                    full.names = TRUE)
mate1.out <- paste(trimmed.reads.dir,(list.files(path = reads.path,
                                                 pattern = read1.pattern,
                                                 full.names = FALSE)), sep = "/")

mate2 <- list.files(path = reads.path,
                    pattern = read2.pattern,
                    full.names = TRUE)
mate2.out <- paste(trimmed.reads.dir,(list.files(path = reads.path,
                                                 pattern = read2.pattern,
                                                 full.names = FALSE)), sep = "/")

# Single end
run_cutadapt(mate1 = mate1,
             mate1.out = mate1.out,
             quality = 25,
             minimum = 17,
             trim.only = TRUE,
             cut = 1,
             adapter1 = "AGATCGGAAGAGCACACGTCT",
             cutadapt = cutadapt.path)

# Paired end
run_cutadapt(mate1 = mate1,
             mate2 = mate2,
             mate1.out = mate1.out,
             mate2.out = mate2.out,
             quality = 25,
             minimum = 17,
             trim.only = TRUE,
             cut = 1,
             adapter1 = "AGATCGGAAGAGCACACGTCT",
             adapter1 = "AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGTAGATCTCGGTGGTCGCCGTATCATT",
             cutadapt = cutadapt.path)
 
## End(Not run)

GrahamHamilton/pipelineTools documentation built on June 19, 2021, 1:08 p.m.