-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathpipeline_mothur_ASVs.txt
109 lines (71 loc) · 4.88 KB
/
pipeline_mothur_ASVs.txt
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
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
#Hacer el stability.files
make.file(inputdir=., type=fastq, prefix=stability)
#hacer contigs
make.contigs(file=stability.files, processors=23)
#Checar los contigs
summary.seqs(fasta=stability.trim.contigs.fasta,count=stability.contigs.count_table)
##### Del SOP
summary.seqs(fasta=stability.trim.contigs.fasta)
make.contigs(file=stability.files, maxambig=0, maxlength=451)
summary.seqs(fasta=stability.trim.contigs.fasta, count=stability.contigs.count_table)
####
#Screen.seqs Recuerda cambiar a tus valores entre el intervalo 2.5 (starts) y 97.5% (end)
screen.seqs(fasta=stability.trim.contigs.fasta, count=stability.contigs.count_table, maxambig=0, maxlength=451)
#mothur, recuerda lo que hacemos
get.current()
#¿qué sobró?
summary.seqs()
#secuencuas unicas
unique.seqs(fasta=stability.trim.contigs.fasta, count=stability.contigs.count_table)
summary.seqs(count=stability.trim.contigs.count_table)
#tabla de cuentas
#count.seqs(name=stability.trim.contigs.good.names, group=stability.contigs.good.groups)
summary.seqs(count=stability.trim.contigs.good.count_table)
###Por alguna razón, Félix se salta el pcr.seqs directo al align.seqs
#Alineamiento de secuencias
align.seqs(fasta=stability.trim.contigs.unique.fasta, reference=silva.full_v138.fasta) # Cambié "silva.v4.fasta" por "silva.full_v138.fasta" completa de 2021 #https://github.com/NorwegianVeterinaryInstitute/BioinfTraining/blob/master/Building_classification_databases.md | y | https://forum.mothur.org/t/fasta-mothur-eft/3409 | y |
#qué salió?
summary.seqs(fasta=stability.trim.contigs.unique.align, count=stability.trim.contigs.count_table)
#Otro screen.seqs Recuerda cambiar a tus valores entre el intervalo 2.5 y 97.5%
screen.seqs(fasta=stability.trim.contigs.good.unique.align, count=stability.trim.contigs.good.count_table, summary=stability.trim.contigs.good.unique.summary, start=6388, end=25316, maxhomop=8)
#qué salió?
summary.seqs(fasta=current, count=current)
#filtrar los gaps y otros espurios
filter.seqs(fasta=stability.trim.contigs.good.unique.good.align, vertical=T, trump=.)
#secuencias unicas
unique.seqs(fasta=stability.trim.contigs.unique.good.filter.fasta, count=stability.trim.contigs.good.count_table)
#precluster
pre.cluster(fasta=stability.trim.contigs.unique.good.filter.unique.fasta, count=stability.trim.contigs.unique.good.filter.count_table, diffs=2)
#chimeras
chimera.vsearch(fasta=stability.trim.contigs.unique.good.filter.unique.precluster.fasta, count=stability.trim.contigs.unique.good.filter.unique.precluster.count_table, dereplicate=t)
#remover chimeras
remove.seqs(fasta=stability.trim.contigs.good.unique.good.filter.unique.precluster.fasta, accnos=stability.trim.contigs.good.unique.good.filter.unique.precluster.denovo.vsearch.accnos)
#qué salió?
summary.seqs(fasta=current, count=current)
#clasificacion taxo
classify.seqs(fasta=stability.trim.contigs.good.unique.good.filter.unique.precluster.pick.fasta, count=stability.trim.contigs.good.unique.good.filter.unique.precluster.denovo.vsearch.pick.count_table, reference=trainset9_032012.pds.fasta, taxonomy=trainset9_032012.pds.tax, cutoff=80)
#remoción de indeseables
remove.lineage(fasta=stability.trim.contigs.good.unique.good.filter.unique.precluster.pick.fasta, count=stability.trim.contigs.good.unique.good.filter.unique.precluster.denovo.vsearch.pick.count_table, taxonomy=stability.trim.contigs.good.unique.good.filter.unique.precluster.pick.pds.wang.taxonomy, taxon=Chloroplast-Mitochondria-unknown-Eukaryota)
#qué salió?
summary.tax(taxonomy=current, count=current)
##Felix se salta los errores del mock
###Pasos para ASV
#make.shared
make.shared(count=stability.trim.contigs.good.unique.good.filter.unique.precluster.denovo.vsearch.pick.pick.count_table)
#Clasifiacipn en ASVs (Esto se salta todos los paso de OTUS)
classify.otu(list=stability.trim.contigs.good.unique.good.filter.unique.precluster.denovo.vsearch.pick.pick.ASV.list, count=stability.trim.contigs.good.unique.good.filter.unique.precluster.denovo.vsearch.pick.pick.count_table, taxonomy=stability.trim.contigs.good.unique.good.filter.unique.precluster.pick.pds.wang.pick.taxonomy, label=ASV)
#renombrar
rename.file(taxonomy=stability.trim.contigs.good.unique.good.filter.unique.precluster.denovo.vsearch.pick.pick.ASV.ASV.cons.taxonomy, shared=stability.trim.contigs.good.unique.good.filter.unique.precluster.denovo.vsearch.pick.pick.ASV.shared)
#contar las secuencias de cada muestra
count.groups(shared=stability.shared)
#Felix se salta el subsample
#rarefaccion
rarefaction.single(stability.shared, calc=sobs, freq=100)
#indices
summary.single(shared=stability.shared, calc=nseqs-coverage-sobs-invsimpson-shannon-chao-simpson-simpsoneven-bootstrap-goodscoverage, subsample=T)
#biom aqui lo deja Félix
make.biom(shared=stability.shared, constaxonomy=stability.taxonomy)
# Beta div
dist.shared(shared=stability.shared, calc=thetayc-jclass, subsample=t)
#pcoa
pcoa(phylip=stability.thetayc.ASV.lt.ave.dist)