BMC

Alias for the type BioMarkovChain.

struct BioMarkovChain{S<:DataType, M<:AbstractMatrix, I<:AbstractVector, N<:Integer} <: AbstractBioMarkovChain

A BioMarkovChain represents a Markov chain used in biological sequence analysis. It contains a transition probability matrix (tpm) and an initial distribution of probabilities (inits) and also the order of the Markov chain.

Fields

• alphabet::A: The state space of the sequence whether DNA, RNA AminoAcid DataTypes.
• tpm::M: The transition probability matrix.
• inits::I: The initial distribution of probabilities.
• n::N: The order of the Markov chain.

Constructors

• BioMarkovChain(tpm::M, inits::I, n::N=1) where {M<:AbstractMatrix, I<:AbstractVector, N<:Integer}: Constructs a BioMarkovChain object with the provided transition probability matrix, initial distribution, and order.
• BioMarkovChain(sequence::LongNucOrView{4}, n::Int64=1): Constructs a BioMarkovChain object based on the DNA sequence and transition order.

Example

sequence = LongDNA{4}("ACTACATCTA")

model = BioMarkovChain(sequence, 2)

BioMarkovChain of DNAAlphabet{4}() and order 2:
  - Transition Probability Matrix -> Matrix{Float64}(4 × 4):
   0.4444  0.1111  0.0     0.4444
   0.4444  0.4444  0.0     0.1111
   0.0     0.0     0.0     0.0
   0.1111  0.4444  0.0     0.4444
  - Initial Probabilities -> Vector{Float64}(4 × 1):
   0.3333  0.3333  0.0     0.3333

initials(sequence::SeqOrView{A}) where A

Calculate the estimated initial probabilities for a Markov chain based on a given sequence.

This function takes a sequence of states and calculates the estimated initial probabilities of each state in the sequence for a Markov chain. The initial probabilities are estimated by counting the occurrences of each state at the beginning of the sequence and normalizing the counts to sum up to 1.

\[\begin{align} \pi{i} &= P(X_{i} = i), i \in T \\ \sum_{i=1}^{N} \pi_{i} &= 1 \end{align}\]

Now using the dinucleotides counts estimating the initials would follow:

\[\hat{\pi_{i}} = c_{i} \sum_{k} c_{k}\]

Arguments

• sequence::SeqOrView{A}: The sequence of states representing the Markov chain.

Returns

An Vector{Flot64} of estimated initial probabilities for each state in the sequence.

log_odds_ratio_matrix(model1::BioMarkovChain, model2::BioMarkovChain)

Calculates the log-odds ratio between the transition probability matrices of two BioMarkovChain models.

\[\beta = \log \frac{P(x|\mathscr{m}_{1})}{P(x|\mathscr{m}_{2})}\]

Where $\mathscr{m}_{1}$ and $\mathscr{m}_{2}$ are the two models transition probability matrices.

Arguments

• model1::BioMarkovChain: The first BioMarkovChain model.
• model2::BioMarkovChain: The second BioMarkovChain model.

log_odds_ratio_score(sequence::SeqOrView{A}, model::BioMarkovChain; b::Number = ℯ)

Compute the log odds ratio score between a given sequence and a BioMarkovChain model.

\[S(x) = \sum_{i=1}^{L} \beta_{x_{i}x} = \sum_{i=1} \log \frac{a^{\mathscr{m}_{1}}_{i-1} x_i}{a^{\mathscr{m}_{2}}_{i-1} x_i}\]

Arguments

• sequence::SeqOrView{A}: A sequence of elements of type A.
• model::BioMarkovChain: A BioMarkovChain model.
• b::Number = ℯ: The base of the logarithm used to compute the log odds ratio.

Returns

The log odds ratio score between the sequence and the model.

Example

markovprobability(sequence::LongNucOrView{4}, model::BioMarkovChain)

Compute the probability of a given sequence using a transition probability matrix and the initial probabilities distributions of a BioMarkovModel.

\[P(X_1 = i_1, \ldots, X_T = i_T) = \pi_{i_1}^{T-1} \prod_{t=1}^{T-1} a_{i_t, i_{t+1}}\]

Arguments

• sequence::LongNucOrView{4}: The input sequence of nucleotides.

Keywords

• model::BioMarkovChain=ECOLICDS: A given BioMarkovChain model.
• logscale::Bool=false: If true, the function will return the log2 of the probability.
• b::Number=2: The base of the logarithm used to compute the log odds ratio.

Returns

• probability::Float64: The probability of the input sequence given the model.

Example

seq = LongDNA{4}("CGCGCGCGCGCGCGCGCGCGCGCGCG")
   
markovprobability(seq, model=CPGPOS, logscale=true)
    -45.073409957110556

markovprobability(seq, model=CPGNEG, logscale=true)
    -74.18912168395339

perronfrobenius(sequence::SeqOrView{A}, n::Int64=1) where A

Compute the Perron-Frobenius matrix, a column-stochastic version of the transition probability matrix (TPM), for a given nucleotide sequence.

The Perron-Frobenius matrix captures the asymptotic probabilities of transitioning between nucleotides in the sequence over a specified number of steps n. It provides insight into the long-term behavior of a Markov chain or a dynamical system associated with the sequence.

Arguments

• sequence::SeqOrView{A}: A nucleotide sequence represented as a NucleicSeqOrView{A} object.
• n::Int64=1: The number of steps to consider for the transition probability matrix. Default is 1.

Returns

A copy of the Perron-Frobenius matrix. Each column of this matrix corresponds to the probabilities of transitioning from the current nucleotide state to all possible nucleotide states after n steps.

Example

sequence = LongSequence{DNAAlphabet{4}}("ACGTCGTCCACTACGACATCAGC")  # Replace with an actual nucleotide sequence
n = 2
pf = perronfrobenius(sequence, n)

transition_count_matrix(sequence::LongSequence{DNAAlphabet{4}})

Compute the transition count matrix (TCM) of a given DNA sequence.

Arguments

• sequence::LongSequence{DNAAlphabet{4}}: a LongSequence{DNAAlphabet{4}} object representing the DNA sequence.

Returns

A Matrix object representing the transition count matrix of the sequence.

Example

seq = LongDNA{4}("AGCTAGCTAGCT")

tcm = transition_count_matrix(seq)

4×4 Matrix{Int64}:
 0  0  3  0
 0  0  0  3
 0  3  0  0
 2  0  0  0

transition_probability_matrix(sequence::LongSequence{DNAAlphabet{4}}, n::Int64=1)

Compute the transition probability matrix (TPM) of a given DNA sequence. Formally it construct $\hat{\mathscr{M}}$ where:

\[\mathscr{m}_{ij} = P(X_t = j \mid X_{t-1} = i) = \frac{{P(X_{t-1} = i, X_t = j)}}{{P(X_{t-1} = i)}}\]

The transition matrices of DNA and Amino-Acids are arranged sorted and in row-wise matrices:

First the DNA matrix:

\[\mathscr{M}_{DNA} = \begin{bmatrix} _{AA} & _{AC} & _{AG} & _{AT} \\ _{CA} & _{CC} & _{CG} & _{CT} \\ _{GA} & _{GC} & _{GG} & _{GT} \\ _{TA} & _{TC} & _{TG} & _{TT} \\ \end{bmatrix}\]

And then, the Aminoacids:

\[\mathscr{M}_{AA} = \begin{bmatrix} _{AA} & _{AC} & _{AD} & \dots & _{AW} \\ _{CA} & _{CC} & _{CD} & \dots & _{CW} \\ _{DA} & _{DC} & _{DD} & \dots & _{DW} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ _{WA} & _{WC} & _{WD} & \dots & _{WW} \\ \end{bmatrix}\]

Arguments

• sequence::LongNucOrView{4}: a LongNucOrView{4} object representing the DNA sequence.
• n::Int64=1: The order of the Markov model. That is the $\hat{M}^{n}$

Keywords

• extended_alphabet::Bool=false: If true will pass the extended alphabet of DNA to search

Returns

A Matrix object representing the transition probability matrix of the sequence.

Example

seq = dna"AGCTAGCTAGCT"

tpm = transition_probability_matrix(seq)

4×4 Matrix{Float64}:
 0.0  0.0  1.0  0.0
 0.0  0.0  0.0  1.0
 0.0  1.0  0.0  0.0
 1.0  0.0  0.0  0.0