| max_mac_ll.h.html | mathcode2html |
| Source file: max_mac_ll.h | |
| Converted: Sun Aug 7 2016 at 13:47:28 | |
| This documentation file will not reflect any later changes in the source file. |
$$\phantom{******** If you see this on the webpage then the
browser could not locate *********}$$
$$\phantom{******** jsMath/easy/load.js or the variable root
is set wrong in this file *********}$$
$$\newcommand{\vector}[1]{\left[\begin{array}{c} #1 \end{array}\right]}$$
$$\newenvironment{matrix}{\left[\begin{array}{cccccccccc}} {\end{array}\right]}$$
$$\newcommand{\A}{{\cal A}}$$
$$\newcommand{\W}{{\cal W}}$$
/* Copyright 2008-2013 Research Foundation State University of New York */
/* This file is part of QUB Express. */
/* QUB Express is free software; you can redistribute it and/or modify */
/* it under the terms of the GNU General Public License as published by */
/* the Free Software Foundation, either version 3 of the License, or */
/* (at your option) any later version. */
/* QUB Express 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 General Public License for more details. */
/* You should have received a copy of the GNU General Public License, */
/* named LICENSE.txt, in the QUB Express program directory. If not, see */
/* <http://www.gnu.org/licenses/>. */
#ifndef MAX_MAC_LL_H
#define MAX_MAC_LL_H
#include "maxmll.h"
#include "../qubfast/qubx_model.h"
#include "callbk_reportfun.h"
#ifdef __cplusplus
extern "C" {
#endif
/*
|
|
Maximum Macroscopic Likelihood[from (Milescu 2005)] by Chris Nicolai 2009This function computes the log-likelihood of data given a model. See also: Up: Index
The ModelOur Markov model is a graph with colored vertices. A vertex is called a "state," and its color is a nonnegative integer called its "class." States in the same class are indistinguishable (same amp). To describe the vertices, provide the array
int clazz[Ns] = [class of each state]
int Ns = number of states
Each edge is labeled with its transition rate (probability per second). These form the matrix
Q, a Ns x Ns matrix with Each \(Q_{a,b} = K0_{a,b} * Ligand_{a,b} * e^{K1_{a,b} * Voltage_{a,b}}\). You provide the Ns x Ns matrices
double **K0, **K1 of kinetic parameters
int **Ligand, **Voltage index of the ligand or voltage signal influencing each rate, or 0
with the diagonals undefined. A pair of states is either connected in both directions or neither. To indicate un-connectedness, set \(K0_{a,b} = K0_{b,a} = 0\) For constant state entry probabilities, as described in (Qin...1996), provide the array int P0[Ns] For equilibrium entry probabilities, as described in (Milescu...2005), provide P0 = NULL TODO: talk about amp/cond model and channel count
The DataThis is wrong -- leftover from MSL -- The data consist of one or more parallel signals. The first (index 0) is the Markovian one to be analyzed. Each additional signal describes a ligand or voltage variable. Signals are in model order; i.e., if (v = Voltage[a][b]) != 0, then the signal at index v is the voltage controlling a->b. For each signal you provide int DwellCount number of events int Classes[DwellCount] class index of each event float Durations[DwellCount] length of each dwell in milliseconds double Amps[ClassCount] amplitude (or ligand/voltage value) of each class and you give the signals together as int Nsignal number of signals int *dwellCounts int **classeses float **durationses double **ampses You provide one or more segments of data, each with the same number of signals. All together they are given as int Nseg int Nsignal int **dwellCountses [segment][signal] int ***classeseses float ***durationseses double ***ampseses We multiplex the signals to create a single idealized signal which changes class whenever any source signal changes class. Each plexi-class denotes a Markov class with a specific set of experimental constants (the other signals' idealized amplitudes). Then, as in MIL, we subtract tdead. tdead is the longest duration of events that can't be reliably detected. MSL deletes any such events, by merging them with their prior. Then for computational reasons, tdead is subtracted from each event, and they're converted to seconds.Sadly, someone has to allocate memory to store this processed signal. It has the form int Nseg int Nsignal int dwellCounts[Nseg] int **classeses float **durationses int Nplex int plexiclasses[2*Nplex] alternating markov-class[i/2], stimclass[i/2] int Nstim double stimclasses[Nstim*Nsignal] [stimcls * Nsignal + signal_ix], signal 0 undefined First, we can scan your data to compute upper bounds for dwellCounts, Nplex and Nstim. You allocate dwellCounts[Nseg] and call: |
*/ // MAXMLL_API void __stdcall mac_multiplex_bounds(int Nseg, int Nsignal, // int **dwellCountses, int ***classeseses, // float ***durationseses, double ***ampseses, // int *out_dwellCounts, int *out_Nplex, int *out_Nstim); // /*
|
|
Then you allocate classeses[Nseg][dwellCounts[i]] and durationses, plexiclasses and stimclasses, and call this function to multiplex and process the signals: |
*/
// MAXMLL_API void __stdcall mac_multiplex(int Nseg, int Nsignal,
// int **dwellCountses, int ***classeseses,
// float ***durationseses, double ***ampseses,
// double tdead_ms, int *out_dwellCounts,
// int **out_classeses, float **out_durationses,
// int *out_Nplex, int *out_plexiclasses,
// int *out_Nstim, double *out_stimclasses);
typedef struct
{
float time, cur, var;
} qubx_mac_sample;
typedef struct
{
float fit, fit_var, ll;
} qubx_mac_sample_out;
// opaque type qubx_mac_data = (void *):
// stores data (uneven sampling possible, or no data) and stimulus idealization
MAXMLL_API void* __stdcall qubx_mac_data_create(int Nstimsig, int data_cap);
MAXMLL_API void __stdcall qubx_mac_data_reset(void *data, int Nstimsig, int data_cap);
MAXMLL_API void __stdcall qubx_mac_data_free(void *data);
MAXMLL_API int __stdcall qubx_mac_data_append_segment(void *data);
MAXMLL_API void __stdcall qubx_mac_data_append_chunk(void *data, int seg_ix, double dur_sec, int Nd,
float *T, float *I, float *V, int *dwellCounts,
int **classeses, float **durationses, double **ampses);
MAXMLL_API double* __stdcall qubx_mac_data_get_stimclasses(void *data, int *Nstimclass);
// returns double[Nstimclass * Nsignal], Nsignal = Nstimsig + 1;
MAXMLL_API double* __stdcall qubx_mac_data_get_stimclass_frac(void *data); // after get_stimclasses
// returns double[Nstimclass]; weight of each stimclass; sums to 1.0; valid after calling get_stimclasses()
MAXMLL_API qubx_mac_sample_out* __stdcall qubx_mac_data_get_buffer(void *data);
typedef struct {
void *storage;
// you set:
void *data;
int useVar;
int optNchannel;
int accel;
// and add_model()s
// ready after setup_constraints():
int Npar, NfPar;
double *pars, *fPars;
// to gray out the GPU checkbox
int gpu_fail;
// output after qubx_mac_ll()
double ll;
// this is what's actually used (initialized from fastmodel nchannel)
double float_Nchannel;
} qubx_mac_param;
MAXMLL_API qubx_mac_param* __stdcall qubx_mac_param_create(int Nmodel);
MAXMLL_API void __stdcall qubx_mac_param_free(qubx_mac_param *p);
MAXMLL_API void __stdcall qubx_mac_param_add_model(qubx_mac_param *p, qubx_model *model, qubx_stim_amps *ampm, qubx_stim_rates *ratesm);
MAXMLL_API int __stdcall qubx_mac_setup_constraints(qubx_mac_param *p);
MAXMLL_API void __stdcall qubx_mac_do_fPars_to_model(qubx_mac_param *p);
MAXMLL_API int __stdcall qubx_mac_ll(qubx_mac_param *p);
MAXMLL_API void __stdcall qubx_mac_calc_std(qubx_mac_param *p, double *InvHessian);
#ifdef __cplusplus
}
#endif
#endif