1. Open the Species selector and enter human in the search box or click on the Human icon in the species list at the bottom of the page. Select to add the reference (GRCh38.p14) assembly and click on the green Add button. You should now see your selected species at the top of the page. Click on Find gene next to the species name, enter MYH9 in the search bar and click Go. In the results, click on MYH9 ENSG00000100345.23 and select Genome Browser. Click on the gene in the genome browser or click on the three dots (…) next to MYH9 protein_coding in the track panel on the right.

   • The gene is located on chromosome 22 on the reverse strand.

   In the track panel, click on +22 transcripts to expand the list of transcripts.

   • MYH9 has 23 transcripts. 6 of these are protein_coding (this includes the MANE select transcript) and 3 are defined as protein_coding_CDS_not_defined (a transcript that belongs to a protein_coding gene and does not contain an open reading frame).

   In the track panel, click on the three dots (…) next to the MANE Select transcript to open the transcript panel and find out more details.

   • The MANE Select is a default transcript per human gene that is representative of biology, well-supported, expressed and highly-conserved. The MANE Select transcript has 40 coding exons.

   In the transcript panel, expand the Download option.

   • You can download the genomic sequence and exons of the gene, and genomic sequence, cDNA (the sequence of the spliced exons of a transcript expressed in DNA notation – T rather than U – representing the coding or sense strand), exons, protein sequence and coding sequence (CDS) of the tanscript.

2. To open the Entity Viewer, scroll to the bottom of the page and click on the Entity Viewer icon. Alternatively, if you are in the genome browser view, click on the first MYH9 transcript and click on View in Entity Viewer in the pop-up menu. Once in the Entity Viewer, click on the MANE Select transcript (this will open a grey information panel underneath).

   • The protein is 1,960 aa long. The Ensembl protein ID is ENSP00000216181.6.

   In the Entity Viewer, switch to the Gene function tab panel, click on the three dots (…) next to the MANE Select transcript to open the transcript panel and find out more details.

   • Protein domains are distinct functional and/or structural units in a protein, which are usually responsible for a particular function or interaction, contributing to the overall role of a protein. There are several methods and algorithms that can be used to classify protein domains into families and functional sites. For the ENSP00000216181.6 protein, PANTHER and Pfam annotations are available.

   In the Gene function tab, click on the UniProtKB/Swiss-Prot P35579 link to open the corresponding entry in UniProt.

   • According to UniProtKB, the function of the protein is as follows: _Cellular myosin that appears to play a role in cytokinesis, cell shape, and specialized functions such as secretion and capping. Required for cortical actin clearance prior to oocyte exocytosis; promotes cell motility in conjunction with S100A4; and during cell spreading, plays an important role in cytoskeleton reorganization, focal contact formation (in the margins but not the central part of spreading cells), and lamellipodial retraction; this function is mechanically antagonized by MYH10.

3. In the Entity Viewer, switch to the Gene relationships tab to view any homologues of the MYH9 gene. You can filter the homology table by header. The % protein similarity is the precentage of identical amino acid residues aligned against each other. Click on the % Protein similarity header once to sort the table in descending order.

   • The Pan troglodytes (Chimpanzee) homologue is the most similar in terms of protein similarity. The gene ID of the homologue is ENSPTRG00000014309.6.

   Click on the gene ID ENSPTRG00000014309.6, to open the corresponding homologue in the Entity Viewer. Count the number of transcripts you can see under the Transcripts tab.

   • The chimp homologue has 4 transcripts.

1. In the Genome Browser app, click on the Navigate browser image icon in the track panel on the right. Expand Go to new location, select Chr 2, enter 176,092,720 in Start and 176,190,907 in End. Click on the green Go button. Count the number of HOXD genes in the genome browser on the left.

   There are 9 HOXD genes on the forward and 1 HOXD gene on the reverse strand in this region.

   Go to the GeneCards entry for the HOXD cluster. You can find diseases associate with the cluster in the Summaries for HOXD@ Gene section.

   • The HOXD gene cluster is associated with synpolydactyly and clubfoot.

2. In the genome browser, click on the HOXD1 (ENSG00000128645.15) gene and select Genome Browser in the pop-up menu. Scroll down to find the Regulatory track, which is indicated with the letter R on the far left. Click on the regulatory feature to see what type it is. Alternatively, you can find a legend in the track panel under the Regulation tab on the right.

   The gene overlaps 1 regulatory feature: the promoter ENSR00000629037.

   In the Regulation (R) track in the genome browser, click on the purple regulatory feature. The pop-up menu will tell you what type of regulatory feature this is. To find a description of the feature, go to the Regulation tab in the track panel on the right. Click on the three dots (…) next to TF binding.

   • This regulatory feature is a TF binding site. These are sites that are enriched for Transcription Factor (TF) binding, but they lack epigenomic evidence to be classified as an enhancer or promoter

3. In the genome browser, scroll down to the Variation (V) track. In the track, you should see green, yellow and salmon coloured variants. You can find a colour legend in the Variation tab in the track panel on the right.

   The HOXD1 region overlaps various transcript, splicing and protein altering variants.

   Using your mouse, zoom in on the first yellow-coloured variant. Click on the first yellow variant to find out more information about it.

   • This variant (rs1691721387) is a single nucleotide variant (SNV). The reference allele is C and the alternative allele is T.

   In the pop-up menu, select View in Entity Viewer to find out more information about the particular variant. In the track panel on the right, you can find the source of the variant (at the top) and the variant location in VCF (at the bottom).

   • This rs1691721387 variant was imported from dbSNP (release 156). The VCF is as follows: 2 176188810 rs1691721387 C T.

1. Open the Species selector and enter bread wheat in the search box or click on the Bread wheat icon in the species list at the bottom of the page. Select to add the IWGSC assembly and click on the green Add button. You should now see your selected species at the top of the page. Click on Find gene next to the species name, enter TraesCS2D02G248400 in the search bar and click Go. In the results, click on TraesCS2D02G248400 and select View in Genome Browser. You can find the number of transcripts in the genome browser view on the left, or the track panel on the right.

   TraesCS2D02G248400 has 2 transcripts.

   Click on the three dots (…) next to the first transcript (TraesCS2D02G248400.2) in the track panel on the right. You can click on the questionmark (?) next to Ensembl canonical to find a description.

   • The Ensembl canonical transcript is a single, representative transcript identified at every locus. The gene codes for an oxygen evolving enhancer protein.

2. Stay in the transcript panel. Expand the Sequences option, select Protein sequence on the right and click on Blast whole sequence. This will take you to the BLAST app. Click on the blue Select species button and select the cultivars Cadenza, Julius and Paragon in Add species. In the BLAST app, click on the green Run button. Click on the blue Results button in the upper right-hand corner to view your results.

   Cadenza has 3, Julius 10 and Paragon 3 hits.

   Click on the Download icon next to the blue Results button in the upper right-hand corner. This will download your results in a compressed folder. Uncompress the folder and open the files using a plain text editor.

   The folder containst 2 files: the output.txt file contains the BLAST alignments, including details about the sequence, alignment scores and statistics, and the sequence alignment. The table.tsv file contains the metadata of your BLAST query, including any BLAST parameters and the results table you saw on the browser.

   Open the output.txt file.

   The Query (top sequence) refers to our sequence of interest (i.e. the protein in the IWGSC genome). The Sbjct (bottom sequence) refers to the homologue (the protein in the cultivar genome).

1. Go to Species selector and search for Escherichia coli by entering the species name in the search bar or by clicking on the species icon in the species list underneath. Select to add the species and click on the green Add button. You should now see the E. coli in your species list at the top of the page. Click on the blue Escherichia coli K-12 ASM584v2 to open the assembly information page. Open the track panel on the right.

   The genome assembly is substrain MG1655.

   Stay in the track panel. You can find the assembly ID under Assembly. Click on the assembly ID to open the corresponding entry in the European Nucleotide Archive (ENA), where the original sequence was submitted to.

   The assembly ID is GCA_000005845.2. In the description in the ENA entry, we can see that MG1655 was derived from strain W1485, which was derived by Joshua Lederberg from the original K-12 isolate obtained from a patient in 1922.

2. In the track panel, click on the Search icon, enter era into the text box and click on era b2566 in the results underneath.

   The era_ gene in the E. coli K-12 genome is also known ‘b2566’ and the coordinates are 2,702,481-2,703,386.

   In the pop-up menu, click on View in Entity Viewer. Switch to the Gene function tab and click on PDBe-KB P06616 to open the corresponding entry in PDBe-KB.

   Accroding to PDBe-KB, the biological function is as follows: An essential GTPase that binds both GDP and GTP, with nucleotide exchange occurring on the order of seconds whereas hydrolysis occurs on the order of minutes. Plays a role in numerous processes, including cell cycle regulation, energy metabolism, as a chaperone for 16S rRNA processing and 30S ribosomal subunit biogenesis. This description is imported from UniProt.

Go to the Ensembl homepage and click on the link Tools at the top of the page. Currently there are nine tools listed in that page. Click on Variant Effect Predictor and enter the three variants as below:

7	117530985	117530985	G/A
7	117531038	117531038	T/C  
7	117531068	117531068	T/C

Note: Variation data input can be done in a variety of formats. See more details about the different data formats and their structure in this VEP documentation page. Click Run. When your job is listed as Done, click View Results.

You will get a table with the consequence terms from the Sequence Ontology project (http://www.sequenceontology.org/) (i.e. synonymous, missense, downstream, intronic, 5’ UTR, 3’ UTR, etc) provided by VEP for the listed SNPs. You can also upload the VEP results as a track and view them on Location pages in Ensembl. SIFT and PolyPhen are available for missense SNPs only. For two of the entered positions, the variations have been predicted to have missense consequences of various pathogenicity (coordinate 117531038 and 117531068), both affecting CFTR. All the three variants have been already annotated and are known as rs1800077, rs1800078 and rs35516286 in dbSNP (databases, literature, etc).

1. Click on VEP at the top of any Ensembl page and open the web interface. Make sure your species is Human. It is good practise to name your VEP jobs something descriptive, such as Patient deletion exercise. Paste the variant in VCF format into the Paste data field and hit Run.

   12 different genes are affected by the SV.

2. Filter your table by select Consequence is transcript_ablation at the top of the table and click Add.

   Yes, there is deletion of complete transcripts of PDS5B, N4BP2L1, BRCA2, RNY1P4, IFIT1P1, ATP8A2P2, N4BP2L2, N4BP2L2-IT2 and one gene without official symbols: ENSG00000212293.

3. To view your variant in the browser click on the location link in the results table 13: 32307062-32908738. The link will open the Region in detail view in a new tab. If you have given your data a name it will appear automatically in red. If not, you may need to Configure this page and add it under the Personal data tab in the pop-up menu.

Exercise 1

You can open the VCF file with any text editor, or by using a command such as: vim 13.32315086-32400268.ALL.chr13_GRCh38.genotypes.20170504.GBR.vcf. Then use the arrow keys to navigate. When done, use the ESC key to switch to normal mode, then write :q! (quit without saving). Lines starting with ## indicate information headers, which contain information about the file content and how it was processed.

Exercise 2

vep -i 13.32315086-32400268.ALL.chr13_GRCh38.genotypes.20170504.GBR.vcf -o VEP_annot_chr13_GBR.vcf --database

Your own script may not look exactly like this and you may employ different flags, here are some common examples:

• --input_file or -i Allows you to specify the location of the input file.
• --output_file or -o Allows you to specify the name of the output file.
• --force_overwrite Allows VEP to overwrite a pre-existing output file with the same name.
• --genomes Points VEP to the Ensembl Genomes (non-vertebrates) server. (not used for this homo sapiens example)
• --database Enable VEP to use local or remote databases.
• --cache Enables the use of the cache (this can speed up VEP significantly - not used here for install simplicity).
• --cache_verson Allows you to specify the cache version. This should be used with Ensembl Genomes caches, since their version numbers do not match Ensembl versions. For example, the VEP/Ensembl version may be 110 and the Ensembl Genomes version 57.
• --check_existing Checks for the existence of known variants that are co-located with your input variants.
• --offline Enables offline mode (no database connections are made).

Exercise 3

You may use a similar less command or method to open the annotated VCF file:

less -S VEP_annot_chr13_GBR.vcf

There are two genes listed in the annotations, ZAR1L ENSG00000189167 and BRCA2 ENSG00000139618. There are multiple rows for each variant representing the variant consequence in different transcripts.

Exercise 4

Rerun VEP with the following command:

vep -i 13.32315086-32400268.ALL.chr13_GRCh38.genotypes.20170504.GBR.vcf -o VEP_annot_chr13_GBR_MANE_Can.txt --canonical --mane --symbol --database

Exercise 5

Filter for BRCA2 and MANE with the following command:

filter_vep -i VEP_annot_chr13_GBR_MANE_Can_gnomadAF.txt -o VEP_annot_FILT.txt --filter "((SYMBOL is BRCA2) and (MANE_SELECT exists))"

By adding transcript set and gene symbol then using filters based on these, we can reduce the annotated variant set down.

To filter to include only the missense results you could add and (Consequence is missense_variant) to the command above.

There is a missense variant - rs80358762 - you can explore phenotype information related to that variant here. Note that although it’s a rare missense variant in BRCA2, the evidence from multiple ClinVar submitters indicate that it is benign.

Bonus!

This variant is novel and does not have an accessioned ID and is not present in the allele frequency resources. Consider just the missense variant in NUDT15, you can see that prediction tools we enabled(alphamissense, revel and CADD) are all relatively low scores, below those that would normally be thought of as pathogenic. However, some ofthe IntAct and MaveDB scores - which are based on experimental data - indicate that the protein structure is significantly impacted, suggesting this variant may merit further investigation.

Exercise 1

You can open the example file with any text editor, or by using a command such as: less example_vars_CHM13T2T.txt. Then use the arrow keys to navigate. When done, use the ESC key to switch to normal mode, then write :q! (quit without saving). There are 5 variants in this example file. rs1230223847 lists two possible allele changes, A > C and A > T.

Exercise 2

vep -i example_vars_CHM13T2T.txt --gff genes.gff3.bgz --fasta softmasked.fa.bgz -o chm13_example_variants.txt --symbol

Your own script may not look exactly like this and you may employ different flags, here are some common examples:

• --input_file or -i Allows you to specify the location of the input file.
• --output_file or -o Allows you to specify the name of the output file.
• --force_overwrite Allows VEP to overwrite a pre-existing output file with the same name.
• --genomes Points VEP to the Ensembl Genomes (non-vertebrates) server. (not used for this homo sapiens example)
• --database Enable VEP to use local or remote databases.
• --cache Enables the use of the cache (this can speed up VEP significantly).
• --cache_verson Allows you to specify the cache version. This should be used with Ensembl Genomes caches, since their version numbers do not match Ensembl versions. For example, the VEP/Ensembl version may be 110 and the Ensembl Genomes version 57.
• --check_existing Checks for the existence of known variants that are co-located with your input variants.
• --offline Enables offline mode (no database connections are made).

Exercise 3

You may use a similar less command or method to open the annotated VCF file:

less -S chm13_example_variants.txt

There are two genes listed in the annotations, ZAR1L ENSG05220054524 and BRCA2 ENSG05220054575. These are different IDs to those in GRCh38. The consequences include: 5_prime_UTR_variant downstream_gene_variant frameshift_variant intron_variant missense_variant splice_acceptor_variant synonymous_variant upstream_gene_variant

Exercise 4

Rerun VEP with the following command:

vep -i example_vars_CHM13T2T.txt --offline --species homo_sapiens_gca009914755v4 --cache_version 107 --fasta chm13_softmasked.fa.bgz --symbol --canonical -o example_vars_CHM13T2T_cache_CAN.txt

Exercise 5

filter_vep -i example_vars_CHM13T2T_cache_CAN.txt --filter "(SYMBOL is BRCA2) and (CANONICAL = YES)" -o example_vars_CHM13T2T_FILT_CAN_BRCA.txt

There are six variants remaining in the file. The multi-allelic rs1230223847 is listed for each of its alleles. Congratulations - you can now annotate variants on the Telomere to Telomere human genome!