mixt 0.1.0

Estimate mixture model weights for a fixed log-likelihood matrix.
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122

// mixt: Riemannian conjugate gradient descent for estimating mixture model weights.
//
// Copyright 2025 Tommi Mäklin [tommi@maklin.fi].
//
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
//
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// Lesser General Public License for more details.
//
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301
// USA
//
use std::path::PathBuf;

use clap::Parser;

// Command-line interface
mod cli;

fn read_counts(
    path: &PathBuf,
    delimiter: u8
) -> Vec<f32> {
    let fs = match std::fs::File::open(path) {
        Ok(fs) => fs,
        Err(e) => panic!("  Error in reading --input-list: {}", e),
    };

    let mut reader = csv::ReaderBuilder::new()
        .delimiter(delimiter)
        .has_headers(false)
        .from_reader(fs);

    reader.records().map(|line| {
        if let Ok(record) = line {
            let val: f32 = (record.iter().next().unwrap().parse::<u32>().unwrap() as f32).ln();
            val
        } else {
            panic!("  Error in reading --weights: {}", path.clone().into_os_string().into_string().unwrap());
        }
    }).collect::<Vec<f32>>()
}

fn read_log_likelihoods(
    path: &PathBuf,
    delimiter: u8
) -> Vec<f32> {
    let fs = match std::fs::File::open(path) {
        Ok(fs) => fs,
        Err(e) => panic!("  Error in reading --input-list: {}", e),
    };

    let mut reader = csv::ReaderBuilder::new()
        .delimiter(delimiter)
        .has_headers(false)
        .from_reader(fs);

    let mut logl: Vec<f32> = Vec::new();
    reader.records().for_each(|line| {
        if let Ok(record) = line {
            record.iter().for_each(|x| { logl.push(x.parse().unwrap()) } );
        } else {
            panic!("  Error in reading --log-likelihood: {}", path.clone().into_os_string().into_string().unwrap());
        }
    });

    logl
}

/// Initializes the logger with verbosity given in `log_max_level`.
fn init_log(log_max_level: usize) {
    stderrlog::new()
        .module(module_path!())
        .quiet(false)
        .verbosity(log_max_level)
        .timestamp(stderrlog::Timestamp::Off)
        .init()
        .unwrap();
}

fn main() {
    let cli = cli::Cli::parse();
    match &cli.command {
        Some(cli::Commands::Fit {
            logl_path,
            weights_path,
            device,
            algorithm,
            tolerance,
            max_iters,
            num_threads,
            verbose,
        }) => {
            init_log(if *verbose { 2 } else { 1 });

            let logl = read_log_likelihoods(logl_path, b'\t');
            let weights = read_counts(weights_path, b'\t');
            let prior: Vec<f32> = vec![1.0; logl.len() / weights.len()];

            let mut options: mixt::OptimizerOpts = Default::default();
            options.tolerance = *tolerance;
            options.max_iters = *max_iters;
            options.device = device.clone().unwrap_or_default();
            options.algorithm = algorithm.clone().unwrap_or_default();

            let (proportions, _) = mixt::optimize_flat(&logl, &weights, &prior, Some(options)).unwrap();

            proportions.iter().enumerate().for_each(|(idx, theta)| {
                eprintln!("{idx}\t{theta}");
            });

        },
        None => {},
    }
}