1 #ifndef INCLUDED_NJN_LOCALMAXSTATUTIL
2 #define INCLUDED_NJN_LOCALMAXSTATUTIL
3 
4 /* $Id: $
5 * ===========================================================================
6 *
7 *                            PUBLIC DOMAIN NOTICE
8 *               National Center for Biotechnology Information
9 *
10 *  This software/database is a "United States Government Work" under the
11 *  terms of the United States Copyright Act.  It was written as part of
12 *  the author's offical duties as a United States Government employee and
13 *  thus cannot be copyrighted.  This software/database is freely available
14 *  to the public for use. The National Library of Medicine and the U.S.
15 *  Government have not placed any restriction on its use or reproduction.
16 *
17 *  Although all reasonable efforts have been taken to ensure the accuracy
18 *  and reliability of the software and data, the NLM and the U.S.
19 *  Government do not and cannot warrant the performance or results that
20 *  may be obtained by using this software or data. The NLM and the U.S.
21 *  Government disclaim all warranties, express or implied, including
22 *  warranties of performance, merchantability or fitness for any particular
23 *  purpose.
24 *
25 *  Please cite the author in any work or product based on this material.
26 *
27 * ===========================================================================*/
28 
29 /*****************************************************************************
30 
31 File name: njn_localmaxstatutil.hpp
32 
33 Author: John Spouge
34 
35 Contents: Random walk parameters
36 
37 ******************************************************************************/
38 
39 #include "njn_matrix.hpp"
40 
41 
42 namespace Njn {
43 	namespace LocalMaxStatUtil {
44 
45 
46         const double REL_TOL = 1.0e-6;
47 
48         void flatten ( // allocates memory for linear probabilities and scores
49         size_t dimension_, // dimension of equilProb_
50         const long int *const *scoreMatrix_, // packed scoring matrix [0...dimension_)[0...dimension_)
51         const double *const *prob_, // prob_ [0...dimension_)[0...dimension_) : distribution of scores sum to 1.0
52         size_t *dim_, // dimension of p_
53         long int **score_, // score [0...dim_) in increasing order
54         double **p_, // linear p_ [0...dim_) : distribution of scores
55         size_t dimension2_ = 0); // dimension2 of equilProb_ : defaults to dimension_
56         // asserts (sum (p_) == 1.0);
57 
58         double lambda (
59         size_t dimMatrix_, // dimension of equilProb_
60         const long int *const *scoreMatrix_, // packed scoring matrix [0...dimension_)[0...dimension_)
61         const double *q_); // q_ [0...dimension_) : distribution of independent letters
62 
63         double mu (
64         size_t dimension_, // #(distinct values)
65         const long int *score_, // scores in increasing order
66         const double *prob_); // probability of corresponding value
67 
68         double lambda (
69         size_t dimension_, // #(distinct values)
70         const long int *score_, // scores in increasing order
71         const double *prob_); // probability of corresponding value
72         // assumes logarithmic regime
73 
74         double muAssoc (
75         size_t dimension_, // #(distinct values) of scores & probabilities (which are paired)
76         const long int *score_, // scores in increasing order
77         const double *prob_, // corresponding probabilities
78         double lambda_ = 0.0); // lambda
79 
80         double muPowerAssoc (
81         size_t dimension_, // #(distinct values) of scores & probabilities (which are paired)
82         const long int *score_, // scores in increasing order
83         const double *prob_, // corresponding probabilities
84         double lambda_ = 0.0, // lambda
85         long int power_ = 1); // power
86 
87         double thetaMin ( // minimizing value for r(theta)
88         size_t dimension_, // #(distinct values)
89         const long int *score_, // scores in increasing order
90         const double *prob_, // probability of corresponding value
91         double lambda_ = 0.0); // lambda
92         // assumes logarithmic regime
93 
94         double rMin ( // minimum value of r(theta)
95         size_t dimension_, // #(distinct values)
96         const long int *score_, // scores in increasing order
97         const double *prob_, // probability of corresponding value
98         double lambda_ = 0.0, // lambda
99         double thetaMin_ = 0.0); // argument of rate
100         // assumes logarithmic regime
101 
102         double r ( // r(theta)
103         size_t dimension_, // #(distinct values)
104         const long int *score_, // scores in increasing order
105         const double *prob_, // probability of corresponding value
106         double theta_); // argument of rate
107         // assumes logarithmic regime
108 
109         long int delta ( // theta [minus delta] for ungapped sequence comparison
110         size_t dimension_, // #(distinct values) of scores & probabilities (which are paired)
111         const long int *score_); // scores
112 
113         double thetaMinusDelta ( // theta [minus delta] for ungapped sequence comparison
114         double lambda_, // lambda, the exponential rate for the local maximum
115         size_t dimension_, // #(distinct values) of scores & probabilities (which are paired)
116         const long int *score_); // scores
117 
118         void descendingLadderEpoch (
119         size_t dimension_, // #(distinct values)
120         const long int *score_, // values
121         const double *prob_, // probability of corresponding value
122         double *eSumAlpha_ = 0, // expectation (sum [alpha])
123         double *eOneMinusExpSumAlpha_ = 0, // expectation [1.0 - exp (sum [alpha])]
124         bool isStrict_ = false, // ? is this a strict descending ladder epoch
125         double lambda0_ = 0.0, // lambda for flattened distribution (avoid recomputation)
126         double mu0_ = 0.0, // mean of flattened distribution (avoid recomputation)
127         double muAssoc0_ = 0.0, // mean of associated flattened distribution (avoid recomputation)
128         double thetaMin0_ = 0.0, // thetaMin of flattened distribution (avoid recomputation)
129         double rMin0_ = 0.0, // rMin of flattened distribution (avoid recomputation)
130         double time_ = 0.0, // get time for the dynamic programming computation
131         bool *terminated_ = 0); // ? Was the dynamic programming computation terminated prematurely ?
132 
133         void descendingLadderEpochRepeat (
134         size_t dimension_, // #(distinct values)
135         const long int *score_, // values
136         const double *prob_, // probability of corresponding value
137         double *eSumAlpha_ = 0, // expectation (sum [alpha])
138         double *eOneMinusExpSumAlpha_ = 0, // expectation [1.0 - exp (sum [alpha])]
139         bool isStrict_ = false, // ? is this a strict descending ladder epoch
140         double lambda_ = 0.0, // lambda for repeats : default is lambda0_ below
141         size_t endW_ = 0, // maximum w plus 1
142         double *pAlphaW_ = 0, // probability {alpha = w} : pAlphaW_ [0, wEnd)
143         double *eOneMinusExpSumAlphaW_ = 0, // expectation [1.0 - exp (sum [alpha]); alpha = w] : eOneMinusExpSumAlphaW_ [0, wEnd)
144         double lambda0_ = 0.0, // lambda for flattened distribution (avoid recomputation)
145         double mu0_ = 0.0, // mean of flattened distribution (avoid recomputation)
146         double muAssoc0_ = 0.0, // mean of associated flattened distribution (avoid recomputation)
147         double thetaMin0_ = 0.0, // thetaMin of flattened distribution (avoid recomputation)
148         double rMin0_ = 0.0, // rMin of flattened distribution (avoid recomputation)
149         double time_ = 0.0, // get time for the dynamic programming computation
150         bool *terminated_ = 0); // ? Was the dynamic programming computation terminated prematurely ?
151         // assumes logarithmic regime
152 
153         bool isProbDist (
154         size_t dimension_, // #(distinct values) of scores & probabilities (which are paired)
155         const double *prob_); // corresponding probabilities
156 
157         bool isScoreIncreasing (
158         size_t dimension_, // #(distinct values)
159         const long int *score_); // scores in increasing order
160 
161         bool isLogarithmic (
162         size_t dimension_, // #(distinct values)
163         const long int *score_, // scores in increasing order
164         const double *prob_); // probability of corresponding value
165 
166 	}
167 }
168 
169 #endif
170 
171