Pairwise Alignment and Dotplot ()
# # Pairwise Alignment and Dotplot # # **Requirements**: # # - Julia 1.6+ # - BioAlignments v2 # - BioSequences v2 # # Pairwise sequence alignment algorithms are provided by the `BioAlignments.jl` package. # The sequence alignment is a kind of optimization problem to find the best scoring alignment. # For example, matching two similar symbols gets a high (often positive) score # and matching dissimilar symbols does a low (often negative) score. using BioSequences using BioAlignments # The following two sequences are homologus protein sequences from [titin](https://en.wikipedia.org/wiki/Titin) # of human and mouse. # Titin is a component of muscle fiber and extremely long (~30,000 amino acids) # compared to other protein sequences. # So, the sequences used here are just part of them. # The complete sequences can be found here: # # - TITIN_HUMAN: <http://www.uniprot.org/uniprot/Q8WZ42> # - TITIN_MOUSE: <http://www.uniprot.org/uniprot/A2ASS6> titin_human = aa""" EYTLLLIEAFPEDAAVYTCEAKNDYGVATTSASLSVEVPEVVSPDQEMPVYPPAIITPLQDTVTSEGQPA RFQCRVSGTDLKVSWYSKDKKIKPSRFFRMTQFEDTYQLEIAEAYPEDEGTYTFVASNAVGQVSSTANLS LEAPESILHERIEQEIEMEMKEFSSSFLSAEEEGLHSAELQLSKINETLELLSESPVYPTKFDSEKEGTG PIFIKEVSNADISMGDVATLSVTVIGIPKPKIQWFFNGVLLTPSADYKFVFDGDDHSLIILFTKLEDEGE YTCMASNDYGKTICSAYLKINSKGEGHKDTETESAVAKSLEKLGGPCPPHFLKELKPIRCAQGLPAIFEY TVVGEPAPTVTWFKENKQLCTSVYYTIIHNPNGSGTFIVNDPQREDSGLYICKAENMLGESTCAAELLVL LEDTDMTDTPCKAKSTPEAPEDFPQTPLKGPAVEALDSEQEIATFVKDTILKAALITEENQQLSYEHIAK ANELSSQLPLGAQELQSILEQDKLTPESTREFLCINGSIHFQPLKEPSPNLQLQIVQSQKTFSKEGILMP EEPETQAVLSDTEKIFPSAMSIEQINSLTVEPLKTLLAEPEGNYPQSSIEPPMHSYLTSVAEEVLSPKEK TVSDTNREQRVTLQKQEAQSALILSQSLAEGHVESLQSPDVMISQVNYEPLVPSEHSCTEGGKILIESAN PLENAGQDSAVRIEEGKSLRFPLALEEKQVLLKEEHSDNVVMPPDQIIESKREPVAIKKVQEVQGRDLLS KESLLSGIPEEQRLNLKIQICRALQAAVA """ #- titin_mouse = aa""" EYTLLLIEAFPEDAAVYTCEAKNDYGVATTSASLSVEVPEVVSPDQEMPVYPPAIVTPLQDTVTSEGRPA RFQCQVSGTDLKVSWYCKDKKIKPSRFFRMTQFEDTYQLEIAEAYPEDEGTYAFVANNAVGQVSSTATLR LEAPESILHERIGQQIEMEMKEIASLLSAEEDFQTYSSDLRLPNANETLELLSEPPARSTQFDSRQEGAA PVFIREISDVEISVEDVAKLSVTVTGCPKPKIQWFFNGMLLTPSADYKFVFDGDTHSLIILFTRFQDEGE YTCLASNEYGKAVCSAHLRISPRGERSTEMESGEKKALEKPKGPCPPYFFKELKPVHCGPGIPAVFEYSV HGEPAPTVLWFKEDMPLYTSVCYTIIHSPDGSGTFIVNDPQRGDSGLYLCKAQNLWGESTCAAELLVLPE DTDVPDASCKEESTLGVPGDFLETSARGPLVQGVDSRQEITAFAEGTISKAALIAEETLQLSYERSVDDS EVGTGVTIGAQKLPPVVLSTPQGTGELPSIDGAVHTQPGRGPPPTLNLQAVQAQTTLPKEATLQFEEPEG VFPGASSAAQVSPVTIKPLITLTAEPKGNYPQSSTAAPDHALLSSVAAETLQLGEKKIPEVDKAQRALLL SQSLAEGCVESLEVPDVAVSNMRSEPQVPFQHTCTEGKILMASADTLKSTGQDVALRTEEGKSLSFPLAL EEKQVLLKEEQSEVVAVPTSQTSKSEKEPEAIKGVKEVREQELLSKETLFPSMPEEQRLHLKTQVRRALQ AAVA """ # The similarity score between two symbols are described by a [substitution matrix](https://en.wikipedia.org/wiki/Substitution_matrix) # and `BioAlignments.jl` provides some predefined matrices commonly seen in bioinformatics. # `BLOSUM62` is a substitution matrix for amino acids # and sometimes used as a default matrix in some programs. # The derivation method and the biological rationale are outlined in this [wikipedia article](https://en.wikipedia.org/wiki/BLOSUM). BLOSUM62 # The score in the matrix can be accessed using a pair of amino acids: BLOSUM62[AA_R, AA_K] # Two protein sequences are aligned to each other using the `pairalign` function, # which takes an alignment optimization problem (`GlobalAlignment()`), # two compared sequences (`titin_human`, `titin_mouse`) and scoring model (`model`). # The algorithms implemented in this module can find the optimal sequence alignments # between two sequences based on affine gap scoring model, # in which inserting gaps of length `k` in a sequence alignment # gets `gap_open + gap_extend * k` where `gap_open` and `gap_extend` # are non-positive scores imposed on starting a new gap and extending a gap, respectively. ## Create a socring model from a substitution matrix and gap scores. model = AffineGapScoreModel(BLOSUM62, gap_open=-11, gap_extend=-1) #- ## Solve a pairwise-alignment optimization problem. result = pairalign(GlobalAlignment(), titin_human, titin_mouse, model) # The achieved score of the alignment is then: score(result) # The global alignment demonstrated above tries to align two sequences in an end-to-end manner. # However, a researcher may be interested in local similarities between two sequences. # In such a case, a [dot plot](https://en.wikipedia.org/wiki/Dot_plot_(bioinformatics)) # is a good visualization to detect local similarities. using CairoMakie """ dotplot(seq1, seq2, submat; windowsize=20, threshold=20) Create a dot plot from `seq1` and `seq2` based on `submat`; `windowsize` specifies the size of a window to accumulate similarity scores along two sequences and `threshold` specifies the score threshold to plot a white dot in a window (when `threshold` is NaN, it plots raw similarity scores). """ function dotplot(seq1, seq2, submat; windowsize=20, threshold=20) m = length(seq1) n = length(seq2) simscores = zeros(m - windowsize + 1, n - windowsize + 1) for i in 1:m-windowsize+1, j in 1:n-windowsize+1 for k in 1:windowsize simscores[i,j] += submat[seq1[i+k-1],seq2[j+k-1]] end end fig = Figure(; resolution=(800,800)) if isnan(threshold) colorscale = :viridis ax = Axis(fig[1,1]; title="dotplot (windowsize=$windowsize, raw similarity score)") else colorscale = :greys ax = Axis(fig[1,1]; title="dotplot (windowsize=$windowsize, threshold=$threshold)") for i in eachindex(simscores) simscores[i] = ifelse(simscores[i] > threshold, 1, 0) end end heatmap!(simscores, colorscale=colorscale) return fig end # We can see a raw similarity score plot by setting `threshold` to `NaN`. dotplot(titin_human, titin_mouse, BLOSUM62, threshold=NaN) # Filtering by the accumulated score highlights similar regions by white dots. # The diagonal breaking lines (from the left-down to the right-up corner) # correspond to the global similarity # between the titin sequences of human and mouse as shown above, # and other shorter lines around it # corerspond to local similarities that are not captured in the global alignment. dotplot(titin_human, titin_mouse, BLOSUM62)
RNA-Seq Coverage (outdated)
# RNA-Seq Coverage # # **Requirements**: # # - Julia 1.6+ # - BioSequences.jl v2 # - BioAlignments.jl v2 # - Xam.jl # - ProgressMeters.jl # # In this tutorial we uses a public RNA-Seq data downloaded from the Sequence Read Archive (SRA): # http://trace.ncbi.nlm.nih.gov/Traces/sra/?run=SRR1238088. # Paired-end FASTQ sequence files were extracted from the archive file # and aligned to a reference genome using a couple of non-julia programs # # ... TODO
The programs under this repository are distributed in the public domain unless otherwise specified.
Source Code: https://github.com/BioJulia/BioTutorials
In this tutorial we’ll learn how to begin programming with R using RStudio. We’ll install R, and RStudio RStudio, an extremely popular development environment for R. We’ll learn the key RStudio features in order to start programming in R on our own.
If you already know how to use RStudio and want to learn some tips, tricks, and shortcuts, check out this Dataquest blog post.
[tidyverse](https://www.dataquest.io/blog/tutorial-getting-started-with-r-and-rstudio/#tve-jump-173bb264c2b)Packages into Memory
#data science tutorials #beginner #r tutorial #r tutorials #rstats #tutorial #tutorials
A famous general is thought to have said, “A good sketch is better than a long speech.” That advice may have come from the battlefield, but it’s applicable in lots of other areas — including data science. “Sketching” out our data by visualizing it using ggplot2 in R is more impactful than simply describing the trends we find.
This is why we visualize data. We visualize data because it’s easier to learn from something that we can see rather than read. And thankfully for data analysts and data scientists who use R, there’s a tidyverse package called ggplot2 that makes data visualization a snap!
In this blog post, we’ll learn how to take some data and produce a visualization using R. To work through it, it’s best if you already have an understanding of R programming syntax, but you don’t need to be an expert or have any prior experience working with ggplot2
#data science tutorials #beginner #ggplot2 #r #r tutorial #r tutorials #rstats #tutorial #tutorials
What exactly is clean data? Clean data is accurate, complete, and in a format that is ready to analyze. Characteristics of clean data include data that are:
Common symptoms of messy data include data that contain:
In this blog post, we will work with five property-sales datasets that are publicly available on the New York City Department of Finance Rolling Sales Data website. We encourage you to download the datasets and follow along! Each file contains one year of real estate sales data for one of New York City’s five boroughs. We will work with the following Microsoft Excel files:
As we work through this blog post, imagine that you are helping a friend launch their home-inspection business in New York City. You offer to help them by analyzing the data to better understand the real-estate market. But you realize that before you can analyze the data in R, you will need to diagnose and clean it first. And before you can diagnose the data, you will need to load it into R!
Benefits of using tidyverse tools are often evident in the data-loading process. In many cases, the tidyverse package
readxl will clean some data for you as Microsoft Excel data is loaded into R. If you are working with CSV data, the tidyverse
readr package function
read_csv() is the function to use (we’ll cover that later).
Let’s look at an example. Here’s how the Excel file for the Brooklyn borough looks:
The Brooklyn Excel file
Now let’s load the Brooklyn dataset into R from an Excel file. We’ll use the
readxlpackage. We specify the function argument
skip = 4 because the row that we want to use as the header (i.e. column names) is actually row 5. We can ignore the first four rows entirely and load the data into R beginning at row 5. Here’s the code:
library(readxl) # Load Excel files brooklyn <- read_excel("rollingsales_brooklyn.xls", skip = 4)
Note we saved this dataset with the variable name
brooklyn for future use.
The tidyverse offers a user-friendly way to view this data with the
glimpse() function that is part of the
tibble package. To use this package, we will need to load it for use in our current session. But rather than loading this package alone, we can load many of the tidyverse packages at one time. If you do not have the tidyverse collection of packages, install it on your machine using the following command in your R or R Studio session:
Once the package is installed, load it to memory:
tidyverse is loaded into memory, take a “glimpse” of the Brooklyn dataset:
glimpse(brooklyn) ## Observations: 20,185 ## Variables: 21 ## $ BOROUGH <chr> "3", "3", "3", "3", "3", "3", "… ## $ NEIGHBORHOOD <chr> "BATH BEACH", "BATH BEACH", "BA… ## $ `BUILDING CLASS CATEGORY` <chr> "01 ONE FAMILY DWELLINGS", "01 … ## $ `TAX CLASS AT PRESENT` <chr> "1", "1", "1", "1", "1", "1", "… ## $ BLOCK <dbl> 6359, 6360, 6364, 6367, 6371, 6… ## $ LOT <dbl> 70, 48, 74, 24, 19, 32, 65, 20,… ## $ `EASE-MENT` <lgl> NA, NA, NA, NA, NA, NA, NA, NA,… ## $ `BUILDING CLASS AT PRESENT` <chr> "S1", "A5", "A5", "A9", "A9", "… ## $ ADDRESS <chr> "8684 15TH AVENUE", "14 BAY 10T… ## $ `APARTMENT NUMBER` <chr> NA, NA, NA, NA, NA, NA, NA, NA,… ## $ `ZIP CODE` <dbl> 11228, 11228, 11214, 11214, 112… ## $ `RESIDENTIAL UNITS` <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 2, 1, 1… ## $ `COMMERCIAL UNITS` <dbl> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0… ## $ `TOTAL UNITS` <dbl> 2, 1, 1, 1, 1, 1, 1, 1, 2, 1, 1… ## $ `LAND SQUARE FEET` <dbl> 1933, 2513, 2492, 1571, 2320, 3… ## $ `GROSS SQUARE FEET` <dbl> 4080, 1428, 972, 1456, 1566, 22… ## $ `YEAR BUILT` <dbl> 1930, 1930, 1950, 1935, 1930, 1… ## $ `TAX CLASS AT TIME OF SALE` <chr> "1", "1", "1", "1", "1", "1", "… ## $ `BUILDING CLASS AT TIME OF SALE` <chr> "S1", "A5", "A5", "A9", "A9", "… ## $ `SALE PRICE` <dbl> 1300000, 849000, 0, 830000, 0, … ## $ `SALE DATE` <dttm> 2020-04-28, 2020-03-18, 2019-0…
glimpse() function provides a user-friendly way to view the column names and data types for all columns, or variables, in the data frame. With this function, we are also able to view the first few observations in the data frame. This data frame has 20,185 observations, or property sales records. And there are 21 variables, or columns.
#data science tutorials #beginner #r #r tutorial #r tutorials #rstats #tidyverse #tutorial #tutorials
In my previous role as a marketing data analyst for a blogging company, one of my most important tasks was to track how blog posts performed.
On the surface, it’s a fairly straightforward goal. With Google Analytics, you can quickly get just about any metric you need for your blog posts, for any date range.
But when it comes to comparing blog post performance, things get a bit trickier.
For example, let’s say we want to compare the performance of the blog posts we published on the Dataquest blog in June (using the month of June as our date range).
But wait… two blog posts with more than 1,000 pageviews were published earlier in the month, And the two with fewer than 500 pageviews were published at the end of the month. That’s hardly a fair comparison!
My first solution to this problem was to look up each post individually, so that I could make an even comparison of how each post performed in their first day, first week, first month, etc.
However, that required a lot of manual copy-and-paste work, which was extremely tedious if I wanted to compare more than a few posts, date ranges, or metrics at a time.
But then, I learned R, and realized that there was a much better way.
In this post, we’ll walk through how it’s done, so you can do my better blog post analysis for yourself!
To complete this tutorial, you’ll need basic knowledge of R syntax and the tidyverse, and access to a Google Analytics account.
Not yet familiar with the basics of R? We can help with that! Our interactive online courses teach you R from scratch, with no prior programming experience required. Sign up and start today!
You’ll also need the
stringr packages installed — which, as a reminder, you can do with the
Finally, you will need a CSV of the blog posts you want to analyze. Here’s what’s in my dataset:
post_url: the page path of the blog post
post_date: the date the post was published (formatted m/d/yy)
category: the blog category the post was published in (optional)
title: the title of the blog post (optional)
Depending on your content management system, there may be a way for you to automate gathering this data — but that’s out of the scope of this tutorial!
For this tutorial, we’ll use a manually-gathered dataset of the past ten Dataquest blog posts.
#data science tutorials #promote #r #r tutorial #r tutorials #rstats #tutorial #tutorials
The Jupyter Notebook is an incredibly powerful tool for interactively developing and presenting data science projects. This article will walk you through how to use Jupyter Notebooks for data science projects and how to set it up on your local machine.
First, though: what is a “notebook”?
A notebook integrates code and its output into a single document that combines visualizations, narrative text, mathematical equations, and other rich media. In other words: it’s a single document where you can run code, display the output, and also add explanations, formulas, charts, and make your work more transparent, understandable, repeatable, and shareable.
Using Notebooks is now a major part of the data science workflow at companies across the globe. If your goal is to work with data, using a Notebook will speed up your workflow and make it easier to communicate and share your results.
Best of all, as part of the open source Project Jupyter, Jupyter Notebooks are completely free. You can download the software on its own, or as part of the Anaconda data science toolkit.
Although it is possible to use many different programming languages in Jupyter Notebooks, this article will focus on Python, as it is the most common use case. (Among R users, R Studio tends to be a more popular choice).
To get the most out of this tutorial you should be familiar with programming — Python and pandas specifically. That said, if you have experience with another language, the Python in this article shouldn’t be too cryptic, and will still help you get Jupyter Notebooks set up locally.
Jupyter Notebooks can also act as a flexible platform for getting to grips with pandas and even Python, as will become apparent in this tutorial.
(In fact, this article was written as a Jupyter Notebook! It’s published here in read-only form, but this is a good example of how versatile notebooks can be. In fact, most of our programming tutorials and even our Python courses were created using Jupyter Notebooks).
#data science tutorials #beginner #jupyter #jupyter notebooks #learn python #pandas #python #tutorial #tutorials