<commit_before>.. _getting_started.rst:

Getting Started
===============

This tutorial provides a comprehensive example of how to use the HGVS
package.  Specifically, we'll:

* parse a transcript (c.) variant in MCL1 obtained from dbSNP
* project that variant to genomic coordinates (as a g. variant)
* project it back on to another transcript in the same gene
* deduce the amino acid change for that variant

We'll use `rs201430561
<http://www.ncbi.nlm.nih.gov/projects/SNP/snp_ref.cgi?rs=201430561>`_ in
MCL1. This gene has several transcripts, and therefore the genomic variant
NC_000001.10:g.150550916G>A has several distinct transcript
representations:

========================  ==========================
transcript (c.)           protein (p.)
========================  ==========================
NM_001197320.1:c.281C>T   NP_001184249.1:p.Ser94Phe
NM_021960.4:c.740C>T      NP_068779.1:p.Ser247Phe
NM_182763.2:c.688+403C>T  *intronic*
========================  ==========================

This variant was chosen because it has data in dbSNP for comparison and
because it has an intronic variant to spice up the example.


Parse the variant
-----------------

To parse variants, we need to create a an instance of the
:py:class:`hgvs.parser.Parser`.  Since building the grammar is
computationally expensive, you should be only one instance and use it for
all parsing operations.

>>> import hgvs.parser
>>> hgvsparser = hgvs.parser.Parser()
>>> var_c1 = hgvsparser.parse_hgvs_variant('NM_001197320.1:c.281C>T')

Parsing a variant results in objects that represent the variant (rather
than, say, Python dictionaries). A SequenceVariant object is comprised of
an accession (``ac``), an HGVS sequence ``type`` (c,g,m,n,r,p), and 0 or
more specific sequence changes (``posedit`` -- a POSition and EDIt).

>>> var_c1
SequenceVariant(ac=NM_001197320.1, type=c, posedit=281C>T)

The ``posedit`` is itself an object:

>>> var_c1.posedit
PosEdit(pos=281, edit=C>T, uncertain=False)

The ``pos`` (position) and ``edit`` attributes are also objects that can
represent intervals and more complex edit operations like indels.  The
``uncertain`` flag enables representation of HGVS uncertainty (typically
with parentheses around the uncertain component).

Finally, "stringifying" a variant regenerates an HGVS variant:

>>> str(var_c1)
'NM_001197320.1:c.281C>T'



Create an HGVSMapper instance
-----------------------------

Mapping variants between genomic (g.), transcript (c.), and protein (p.)
sequences is performed by an instance of :py:class:`HGVSMapper`. As with
the parser, you need only one instance per session.

Variant mapping and validation requires access to external data,
specifically exon structures, transcript alignments, and protein
accessions.  Right now, the only source of this data is via the BDI and
UTA sister projects.  (If you want more information on the architecture of
HGVS, BDI, and UTA, see :doc:`architecture.rst`.  However, you don't
really need to understand the architecture to use HGVS.)

First, connect to UTA via BDI:

>>> import bdi.sources.uta0_pg
bdi = bdi.sources.uta0_pg.UTA0()

Then, with that connection, instantiate an HGVSMapper:

>>> import hgvs.hgvsmapper
hgvsmapper = hgvs.hgvsmapper.HGVSMapper(bdi)

We can use this mapper to transform our transcript variant to a protein variant:

>>> hgvsmapper.hgvsc_to_hgvsp(var_c1,None)
SequenceVariant(ac=NP_001184249.1, type=p, posedit=Ser94Phe)


Map our variant to the genome
-----------------------------
<commit_msg>doc updates<commit_after>.. _getting_started.rst:

Getting Started
===============

This tutorial provides a comprehensive example of how to use the HGVS
package.  Specifically, we'll:

* install hgvs
* parse a transcript (c.) variant in MCL1 obtained from dbSNP
* project that variant to genomic coordinates (as a g. variant)
* project it back on to another transcript in the same gene
* deduce the amino acid change for that variant

We'll use `rs201430561
<http://www.ncbi.nlm.nih.gov/projects/SNP/snp_ref.cgi?rs=201430561>`_ in
MCL1. This gene has several transcripts, and therefore the genomic variant
NC_000001.10:g.150550916G>A has several distinct transcript
representations:

========================  ==========================
transcript (c.)           protein (p.)
========================  ==========================
NM_001197320.1:c.281C>T   NP_001184249.1:p.Ser94Phe
NM_021960.4:c.740C>T      NP_068779.1:p.Ser247Phe
NM_182763.2:c.688+403C>T  *intronic*
========================  ==========================

This variant was chosen because it has data in dbSNP for comparison and
because it has an intronic variant to spice up the example.


Installation
------------

For this demo, you'll need hgvs (of course).  We recommend that you
install IPython as well.  In a reasonably modern environment, the
following should suffice::

  $ pip install hgvs ipython

If you get a warning about setuptools versions, do this::

  $ pip install --upgrade setuptools



Parse the variant
-----------------

To parse variants, we need to create a an instance of the
:py:class:`hgvs.parser.Parser`.  Since building the grammar is
computationally expensive, you should be only one instance and use it for
all parsing operations.  Start ``ipython``, then do this:

>>> import hgvs.parser
>>> hgvsparser = hgvs.parser.Parser()
>>> var_c1 = hgvsparser.parse_hgvs_variant('NM_001197320.1:c.281C>T')

Parsing a variant results in objects that represent the variant (rather
than, say, Python dictionaries). A SequenceVariant object is comprised of
an accession (``ac``), an HGVS sequence ``type`` (c,g,m,n,r,p), and 0 or
more specific sequence changes (``posedit`` -- a POSition and EDIt).

>>> var_c1
SequenceVariant(ac=NM_001197320.1, type=c, posedit=281C>T)

The ``posedit`` is itself an object:

>>> var_c1.posedit
PosEdit(pos=281, edit=C>T, uncertain=False)

The ``pos`` (position) and ``edit`` attributes are also objects that can
represent intervals and more complex edit operations like indels.  The
``uncertain`` flag enables representation of HGVS uncertainty (typically
with parentheses around the uncertain component).

Finally, "stringifying" a variant regenerates an HGVS variant:

>>> str(var_c1)
'NM_001197320.1:c.281C>T'



Create an HGVSMapper instance
-----------------------------

Mapping variants between genomic (g.), transcript (c.), and protein (p.)
sequences is performed by an instance of :py:class:`HGVSMapper`. As with
the parser, you need only one instance per session.

Variant mapping and validation requires access to external data,
specifically exon structures, transcript alignments, and protein
accessions.  Right now, the only source of this data is via the BDI and
UTA sister projects.  (If you want more information on the architecture of
HGVS, BDI, and UTA, see :doc:`architecture.rst`.  However, you don't
really need to understand the architecture to use HGVS.)

First, connect to UTA via BDI:

>>> import bdi.sources.uta0
bdi = bdi.sources.uta0.connect()

By default, you'll connect to the public UTA database instance hosted by
`Invitae <http://invitae.com/>`_.

Then, with that connection, instantiate an HGVSMapper:

>>> import hgvs.hgvsmapper
hgvsmapper = hgvs.hgvsmapper.HGVSMapper(bdi)

We can use this mapper to transform our transcript variant to a protein variant:

>>> hgvsmapper.hgvsc_to_hgvsp(var_c1,None)
SequenceVariant(ac=NP_001184249.1, type=p, posedit=Ser94Phe)


Map our variant to the genome
-----------------------------

{'fraction_non_alphanumeric': 0.1464619492656876, 'fraction_numerical': 0.04205607476635514, 'mean_word_length': 4.512141280353201, 'pattern_counts': {'":': 0, '<': 6, '<?xml version=': 0, '>': 76, 'https://': 0, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '3602981', 'n_tokens_mistral': 2555, 'n_tokens_neox': 2183, 'n_words': 854}