jmschrei / pomegranate

Fast, flexible and easy to use probabilistic modelling in Python.
http://pomegranate.readthedocs.org/en/latest/
MIT License
3.29k stars 591 forks source link

Empty viterbi_path #979

Closed malonzm1 closed 1 year ago

malonzm1 commented 2 years ago

Hi!

viterbi_likelihood, viterbi_path = model.viterbi(sequence) returns -inf for viterbi_likelihood and an empty viterbi_path. Below is the code. For a different dataset it worked. total_reads_Y.txt.tar.gz methylated_reads_Y.txt.tar.gz

import numpy as np
from pomegranate import *
import os

def hidden_states(bsTot_list,bsMeth_list):
    ratio = np.divide(bsMeth_list,bsTot_list)
    control = np.nanmean(ratio[:,:2],axis=1)
    case = np.nanmean(ratio[:,2:4],axis=1)
    diff = case - control
    return diff

chroms = ['Y']
sequences = []
for chrom in chroms:
    bsTot_file = 'total_reads_%s.txt'%(chrom)
    bsTot_list = np.loadtxt(bsTot_file,delimiter='\t',skiprows=0,dtype='float', usecols=range(1, 5))
    bsTot_list[bsTot_list < 5] = np.nan
    bsMeth_file = 'methylated_reads_%s.txt'%(chrom)
    bsMeth_list = np.loadtxt(bsMeth_file,delimiter='\t',skiprows=0,dtype='float', usecols=range(1, 5))
    sequence=hidden_states(bsTot_list,bsMeth_list)
    sequences.append(sequence)
s0 = State(NormalDistribution(0, 0.08), name = 's0')
s1 = State(NormalDistribution(0.3, 0.06), name = 's1')
s2 = State(NormalDistribution(-0.3, 0.06), name = 's2')
model = HiddenMarkovModel()
model.add_states(s0, s1, s2)
model.add_transition(model.start, s0, 0.333)
model.add_transition(model.start, s1, 0.333)
model.add_transition(model.start, s2, 0.333)
model.add_transition(s0, s0, 0.5)
model.add_transition(s0, s1, 0.25)
model.add_transition(s0, s2, 0.25)
model.add_transition(s1, s0, 0.25)
model.add_transition(s1, s1, 0.5)
model.add_transition(s1, s2, 0.25)
model.add_transition(s2, s0, 0.25)
model.add_transition(s2, s1, 0.25)
model.add_transition(s2, s2, 0.5)
model.add_transition(s0, model.end, 0.333)
model.add_transition(s1, model.end, 0.333)
model.add_transition(s2, model.end, 0.333)
model.bake()
model.fit(sequences, distribution_inertia=1.0)
viterbi_likelihood, viterbi_path = model.viterbi(sequences[0])

Pls. advise.

jmschrei commented 1 year ago

Thank you for opening an issue. pomegranate has recently been rewritten from the ground up to use PyTorch instead of Cython (v1.0.0), and so all issues are being closed as they are likely out of date. Please re-open or start a new issue if a related issue is still present in the new codebase.