File: loading_sequences.rst

package info (click to toggle)
python-cogent 2024.5.7a1%2Bdfsg-3
  • links: PTS, VCS
  • area: main
  • in suites: sid
  • size: 74,600 kB
  • sloc: python: 92,479; makefile: 117; sh: 16
file content (233 lines) | stat: -rwxr-xr-x 7,442 bytes parent folder | download
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
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
.. jupyter-execute::
    :hide-code:

    import set_working_directory

.. _load_seq:

Loading a sequence from a file
------------------------------

It's also possible to load a sequence from a :ref:`url <load_url>`.

.. jupyter-execute::

    from cogent3 import load_seq

    seq = load_seq("data/mycoplasma-genitalium.fa", moltype="dna")
    seq

.. warning:: If a file has more than one sequence, only the first one is loaded.

.. jupyter-execute::

    seq = load_seq("data/brca1-bats.fasta", moltype="dna")
    seq

.. note:: The filename suffix is used to infer the data format.

.. _load-seqs:

Loading an alignment from a file or url
---------------------------------------

.. author, Gavin Huttley, Tony Walters, Tom Elliott

Loading aligned sequences
^^^^^^^^^^^^^^^^^^^^^^^^^

.. jupyter-execute::

    from cogent3 import load_aligned_seqs

    aln = load_aligned_seqs("data/long_testseqs.fasta", moltype="dna")
    type(aln)

The load functions record the origin of the data in the ``info`` attribute under a `"source"` key.

.. jupyter-execute::

    aln.info.source

.. note:: The function ``load_aligned_seqs()`` returns an ``ArrayAlignment`` by default. If you set the argument ``array_align=False``, you will get an ``Alignment``. (That class can be annotated.)

.. todo:: add cross ref for description of Info class

Loading unaligned sequences
^^^^^^^^^^^^^^^^^^^^^^^^^^^

The ``load_unaligned_seqs()`` function returns a sequence collection.

.. jupyter-execute::

    from cogent3 import load_unaligned_seqs

    seqs = load_unaligned_seqs("data/long_testseqs.fasta", moltype="dna")
    type(seqs)

.. _load_url:

Loading from a url
^^^^^^^^^^^^^^^^^^

The ``cogent3`` load functions support loading from a url. We load the above fasta file directly from GitHub.

.. jupyter-execute::

    from cogent3 import load_aligned_seqs

    aln = load_aligned_seqs(
        "https://raw.githubusercontent.com/cogent3/cogent3/develop/doc/data/long_testseqs.fasta",
        moltype="dna",
    )

Specifying the file format
^^^^^^^^^^^^^^^^^^^^^^^^^^

The loading functions use the filename suffix to infer the file format. This can be overridden using the ``format`` argument.

.. jupyter-execute::

    from cogent3 import load_aligned_seqs

    aln = load_aligned_seqs("data/long_testseqs.fasta", moltype="dna", format="fasta")
    aln

Specifying the sequence molecular type
--------------------------------------

Simple case of loading a ``list`` of aligned amino acid sequences in FASTA format, with and without ``moltype`` specification. When ``moltype`` is not specified it defaults to ``BYTES`` for the ``ArrayAlignment`` class, ``ASCII`` for the ``Alignment`` class.

.. jupyter-execute::

    from cogent3 import make_aligned_seqs

    protein_seqs = {"seq1": "DEKQL-RG", "seq2": "DDK--SRG"}
    proteins_loaded = make_aligned_seqs(protein_seqs)
    proteins_loaded.moltype
    proteins_loaded

.. jupyter-execute::

    proteins_loaded = make_aligned_seqs(protein_seqs, moltype="protein")
    proteins_loaded

.. note:: This applies to both the ``load_*`` or ``make_*`` functions.

Making an alignment from standard python objects
------------------------------------------------

From a dict of strings
^^^^^^^^^^^^^^^^^^^^^^

.. jupyter-execute::

    from cogent3 import make_aligned_seqs

    seqs = {"seq1": "AATCG-A", "seq2": "AATCGGA"}
    seqs_loaded = make_aligned_seqs(seqs)

From a series of strings
^^^^^^^^^^^^^^^^^^^^^^^^

.. jupyter-execute::

    from cogent3 import make_aligned_seqs

    seqs = {"seq1": "AATCG-A", "seq2": "AATCGGA"}
    seqs_loaded = make_aligned_seqs(seqs)
    seqs_loaded

Stripping label characters on loading
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

Load a list of aligned nucleotide sequences, while specifying the DNA molecule type and stripping the comments from the label. In this example, stripping is accomplished by passing a function that removes everything after the first whitespace to the ``label_to_name`` parameter.

.. jupyter-execute::

    from cogent3 import make_aligned_seqs

    DNA_seqs = {
        "sample1 Mus musculus": "AACCTGC--C",
        "sample2 Gallus gallus": "AAC-TGCAAC",
    }
    loaded_seqs = make_aligned_seqs(
        DNA_seqs, moltype="dna", label_to_name=lambda x: x.split()[0]
    )
    loaded_seqs

Making a sequence collection from standard python objects
---------------------------------------------------------

This is done using ``make_unaligned_seqs()``, which returns a ``SequenceCollection`` instance. The function arguments match those of ``make_aligned_seqs()``. We demonstrate only for the case where the input data is a ``dict``.

.. jupyter-execute::

    from cogent3 import make_unaligned_seqs

    seqs = {"seq1": "AATCA", "seq2": "AATCGGA"}
    seqs = make_unaligned_seqs(data=seqs, moltype="dna")
    seqs

Loading sequences using format parsers
--------------------------------------

``load_aligned_seqs()`` and ``load_unaligned_seqs()`` are just convenience interfaces to format parsers. It can sometimes be more effective to use the parsers directly, say when you don't want to load everything into memory.

Loading FASTA sequences from an open file or list of lines
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

To load FASTA formatted sequences directly, you can use the ``MinimalFastaParser``.

.. note:: This returns the sequences as strings.

.. jupyter-execute::

    from cogent3.parse.fasta import MinimalFastaParser

    f = open("data/long_testseqs.fasta")
    seqs = [(name, seq) for name, seq in MinimalFastaParser(f)]
    seqs

Handling overloaded FASTA sequence labels
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

The FASTA label field is frequently overloaded, with different information fields present in the field and separated by some delimiter. This can be flexibly addressed using the ``LabelParser``. By creating a custom label parser, we can decide which part we use as the sequence name. We show how to convert a field into something specific.

.. jupyter-execute::

    from cogent3.parse.fasta import LabelParser

    def latin_to_common(latin):
        return {"Homo sapiens": "human", "Pan troglodtyes": "chimp"}[latin]

    label_parser = LabelParser(
        "%(species)s", [[1, "species", latin_to_common]], split_with=":"
    )
    for label in ">abcd:Homo sapiens:misc", ">abcd:Pan troglodtyes:misc":
        label = label_parser(label)
        print(label, type(label))

``RichLabel`` objects have an ``Info`` object as an attribute, allowing specific reference to all the specified label fields.

.. jupyter-execute::

    from cogent3.parse.fasta import LabelParser, MinimalFastaParser

    fasta_data = [
        ">gi|10047090|ref|NP_055147.1| small muscle protein, X-linked [Homo sapiens]",
        "MNMSKQPVSNVRAIQANINIPMGAFRPGAGQPPRRKECTPEVEEGVPPTSDEEKKPIPGAKKLPGPAVNL",
        "SEIQNIKSELKYVPKAEQ",
        ">gi|10047092|ref|NP_037391.1| neuronal protein [Homo sapiens]",
        "MANRGPSYGLSREVQEKIEQKYDADLENKLVDWIILQCAEDIEHPPPGRAHFQKWLMDGTVLCKLINSLY",
        "PPGQEPIPKISESKMAFKQMEQISQFLKAAETYGVRTTDIFQTVDLWEGKDMAAVQRTLMALGSVAVTKD",
    ]
    label_to_name = LabelParser(
        "%(ref)s",
        [[1, "gi", str], [3, "ref", str], [4, "description", str]],
        split_with="|",
    )
    for name, seq in MinimalFastaParser(fasta_data, label_to_name=label_to_name):
        print(name)
        print(name.info.gi)
        print(name.info.description)