Comments (10)
matsen
commented on June 23, 2024
The master challenge
Say we have an HMM laid out like
*-*-*
| |
*-*-*
emitting binary symbols.
Say the emission probabilities for the top three states (row 0 corresponds to emitting symbol 0, etc) are
0.1, 0.2, 0.3
0.9, 0.8, 0.7
and those on the bottom are
0.11, 0.13, 0.17
0.89, 0.87, 0.83
and the probability of continuing to the upper right state is 0.23.
Q: what are the marginal and Viterbi probabilities of the sequence 0, 1, 1, 0, 0?
from linearham.
matsen
commented on June 23, 2024
I get (note these are "raw" probabilities, not log probs):
Ham test 1 marginal: 0.00149626
Ham test 1 viterbi: 0.00121161
from linearham.
psathyrella
commented on June 23, 2024
maybe this (row 0 corresponds to state 0, etc) --> (column 0 corresponds to state 0, etc)?
from linearham.
psathyrella
commented on June 23, 2024
I get 0.00121 for the viterbi; I don't know what you mean by the marginal probability of the sequence. Does that mean the probability of the sequence, i.e. forward?
from linearham.
psathyrella
commented on June 23, 2024
oh, I guess so, my fwd prob is 0.001496. full output:
viterbi path (log prob -6.7158):
sequence: z o o z z
path: t t b b b
forward log prob: -6.50479
nice!
from linearham.
psathyrella
commented on June 23, 2024
yaml, for future reference:
name: binample
tracks:
bits: [z,o]
states:
- name: init
transitions:
top_left: 1.0
- name: top_left
emissions:
probs:
z: 0.1
o: 0.9
transitions:
top_middle: 1.
- name: top_middle
emissions:
probs:
z: 0.2
o: 0.8
transitions:
top_right: 0.23
bottom_left: 0.77
- name: top_right
emissions:
probs:
z: 0.3
o: 0.7
transitions:
bottom_middle: 1.0
- name: bottom_left
emissions:
probs:
z: 0.11
o: 0.89
transitions:
bottom_middle: 1.0
- name: bottom_middle
emissions:
probs:
z: 0.13
o: 0.87
transitions:
bottom_right: 1.0
- name: bottom_right
emissions:
probs:
z: 0.17
o: 0.83
transitions:
end: 1.0
from linearham.
matsen
commented on June 23, 2024
Hawt dawg!
We'll talk IRL about next steps.
from linearham.
matsen
commented on June 23, 2024
This one's done; more tests to come.
from linearham.
matsen
commented on June 23, 2024
The Ham test, for posterity's sake:
TEST_CASE("Ham Comparison 1", "[ham]") {
Eigen::VectorXd landing_a(3);
landing_a << 1, 1, 1;
Eigen::MatrixXd emission_matrix_a(2,3);
emission_matrix_a <<
0.1, 0.2, 0.3,
0.9, 0.8, 0.7;
Eigen::VectorXd next_transition_a(2);
next_transition_a << 1, 0.23;
Germline germline_a(landing_a, emission_matrix_a, next_transition_a);
Eigen::VectorXi emission_indices_a(3);
emission_indices_a << 0, 1, 1;
Smooshable s_a = SmooshableGermline(germline_a, 0, emission_indices_a, 1, 2);
Eigen::MatrixXd correct_marginal_a(1, 2);
correct_marginal_a << 0.1*0.8*0.77, 0.1*0.8*0.23*0.7;
REQUIRE(s_a.marginal().isApprox(correct_marginal_a));
REQUIRE(s_a.scaler_count() == 0);
Eigen::VectorXd landing_b(3);
landing_b << 1, 1, 1;
Eigen::MatrixXd emission_matrix_b(2,3);
emission_matrix_b <<
0.11, 0.13, 0.17,
0.89, 0.87, 0.83;
Eigen::VectorXd next_transition_b(2);
next_transition_b << 1, 1;
Germline germline_b(landing_b, emission_matrix_b, next_transition_b);
Eigen::VectorXi emission_indices_b(3);
emission_indices_b << 1, 0, 0;
Smooshable s_b = SmooshableGermline(germline_b, 0, emission_indices_b, 2, 1);
Eigen::MatrixXd correct_marginal_b(2, 1);
correct_marginal_b <<
0.89*0.13*0.17,
0.13*0.17;
REQUIRE(s_b.marginal().isApprox(correct_marginal_b));
REQUIRE(s_b.scaler_count() == 0);
Smooshable s_ab;
Eigen::MatrixXi viterbi_idx_ab;
std::tie(s_ab, viterbi_idx_ab) = Smoosh(s_a, s_b);
Eigen::MatrixXd correct_marginal_ab(1, 1);
correct_marginal_ab <<
(0.1*0.8*0.77*0.89*0.13*0.17 + 0.1*0.8*0.23*0.7*0.13*0.17);
REQUIRE(s_ab.marginal().isApprox(correct_marginal_ab));
Eigen::MatrixXd correct_viterbi_ab(1, 1);
correct_viterbi_ab << 0.1*0.8*0.77*0.89*0.13*0.17;
REQUIRE(s_ab.viterbi().isApprox(correct_viterbi_ab));
REQUIRE(s_ab.scaler_count() == 0);
// now let's test for underflow
landing_b.array() *= SCALE_THRESHOLD;
germline_b = Germline(landing_b, emission_matrix_b, next_transition_b);
s_b = SmooshableGermline(germline_b, 0, emission_indices_b, 2, 1);
REQUIRE(s_b.marginal().isApprox(correct_marginal_b));
REQUIRE(s_b.scaler_count() == 1);
std::cout << "Ham test 1 marginal: " << correct_marginal_ab << std::endl;
std::cout << "Ham test 1 viterbi: " << correct_viterbi_ab << std::endl;
}from linearham.
matsen
commented on June 23, 2024
This is now deprecated in favor of more sophisticated bcrham tests.
from linearham.
Related Issues (20)
- remove dockerhub image HOT 1
- git describe crashes build HOT 1
- Do we really want to use scons? HOT 3
- NodeList object where it expected a string HOT 1
- change actions to positional argument
- make default tag on quay.io HOT 1
- Use new partis fcn add_seqs_to_line() HOT 4
- Cause of huge number of low-probability nodes HOT 6
- update readme with newer (?) output files HOT 2
- update partis + maybe add --min-cluster-size arg HOT 1
- Docker container doesn't build HOT 2
- light chain crash HOT 3
- Github Actions & Quay Container building HOT 2
- Update Master -> Main
- crash in asr/ess script HOT 7
- Inferred ancestral sequences have mutations in ambiguous regions
- Branch lengths almost all much larger than hamming distance
- Add option for less verbose output dir HOT 1
- Failed Quay build from trigger
- How to handle paired h/l data
Recommend Projects
-
React
A declarative, efficient, and flexible JavaScript library for building user interfaces.
-
Vue.js
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
-
Typescript
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
-
TensorFlow
An Open Source Machine Learning Framework for Everyone
-
Django
The Web framework for perfectionists with deadlines.
-
Laravel
A PHP framework for web artisans
-
D3
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
-
Recommend Topics
-
javascript
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
-
web
Some thing interesting about web. New door for the world.
-
server
A server is a program made to process requests and deliver data to clients.
-
Machine learning
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
-
Visualization
Some thing interesting about visualization, use data art
-
Game
Some thing interesting about game, make everyone happy.
Recommend Org
-
Facebook
We are working to build community through open source technology. NB: members must have two-factor auth.
-
Microsoft
Open source projects and samples from Microsoft.
-
Google
Google ❤️ Open Source for everyone.
-
Alibaba
Alibaba Open Source for everyone
-
D3
Data-Driven Documents codes.
-
Tencent
China tencent open source team.
from linearham.