source: branches/pars/DIST/di_matr.hxx @ 12957

Last change on this file since 12957 was 12957, checked in by westram, 11 years ago
  • adds test for distance matrix
    • calculates matrices for
      • all sequence types
      • using 3 different transformations
        • autodetected for dna and pro
        • DI_TRANSFORMATION_FELSENSTEIN for rna (triggers overall frequency calculation)
    • save matrices in all available formats
    • covers all DIST-code-locations marked in [12942] (calculation using user-defined transformation-matrix remains untested)
  • few minor fixes
    • param constness
    • free/delete mismatch
    • useless access behind EOS (while printing NDS generated string in DI_SAVE_READABLE-mode)
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 7.0 KB
Line 
1// =============================================================== //
2//                                                                 //
3//   File      : di_matr.hxx                                       //
4//   Purpose   :                                                   //
5//                                                                 //
6//   Institute of Microbiology (Technical University Munich)       //
7//   http://www.arb-home.de/                                       //
8//                                                                 //
9// =============================================================== //
10
11#ifndef DI_MATR_HXX
12#define DI_MATR_HXX
13
14#ifndef AP_PRO_A_NUCS_HXX
15#include <AP_pro_a_nucs.hxx>
16#endif
17#ifndef AP_TREE_HXX
18#include <AP_Tree.hxx>
19#endif
20#ifndef AP_MATRIX_HXX
21#include <AP_matrix.hxx>
22#endif
23#ifndef _GLIBCXX_STRING
24#include <string>
25#endif
26
27#define di_assert(cond) arb_assert(cond)
28
29#define AWAR_DIST_PREFIX           "dist/"
30#define AWAR_DIST_CORR_TRANS       AWAR_DIST_PREFIX "correction/trans"
31#define AWAR_DIST_SAVE_MATRIX_BASE "tmp/" AWAR_DIST_PREFIX "save_matrix"
32
33#define AWAR_DIST_DIST_PREFIX AWAR_DIST_PREFIX "dist/"
34#define AWAR_DIST_MIN_DIST    AWAR_DIST_DIST_PREFIX "lower"
35#define AWAR_DIST_MAX_DIST    AWAR_DIST_DIST_PREFIX "upper"
36
37enum DI_TRANSFORMATION {
38    DI_TRANSFORMATION_NONE,
39    DI_TRANSFORMATION_SIMILARITY,
40    DI_TRANSFORMATION_JUKES_CANTOR,
41    DI_TRANSFORMATION_FELSENSTEIN,
42
43    DI_TRANSFORMATION_PAM,
44    DI_TRANSFORMATION_CATEGORIES_HALL,
45    DI_TRANSFORMATION_CATEGORIES_BARKER,
46    DI_TRANSFORMATION_CATEGORIES_CHEMICAL,
47
48    DI_TRANSFORMATION_KIMURA,
49    DI_TRANSFORMATION_OLSEN,
50    DI_TRANSFORMATION_FELSENSTEIN_VOIGT,
51    DI_TRANSFORMATION_OLSEN_VOIGT,
52    DI_TRANSFORMATION_ML,
53
54    DI_TRANSFORMATION_FROM_TREE,
55
56    // -------------------- real transformations are above
57
58    DI_TRANSFORMATION_COUNT,         // amount of real transformations
59    DI_TRANSFORMATION_NONE_DETECTED, // nothing was auto-detected
60};
61
62enum DI_MATRIX_TYPE {
63    DI_MATRIX_FULL,
64    DI_MATRIX_COMPRESSED
65};
66
67class AP_sequence_parsimony;
68class AP_sequence_simple_protein;
69class DI_MATRIX;
70
71class DI_ENTRY : virtual Noncopyable {
72    DI_MATRIX *phmatrix;
73    char      *full_name;
74
75public:
76    DI_ENTRY(GBDATA *gbd, DI_MATRIX *phmatri);
77    DI_ENTRY(char *namei, DI_MATRIX *phmatri);
78    ~DI_ENTRY();
79
80    // @@@ make a union?
81    AP_sequence                *sequence;
82    AP_sequence_parsimony      *sequence_parsimony; // if exist ok
83    AP_sequence_simple_protein *sequence_protein;
84
85    char *name;
86    int   group_nr;                                 // species belongs to group number xxxx
87};
88
89enum DI_SAVE_TYPE {
90    DI_SAVE_PHYLIP_COMP,
91    DI_SAVE_READABLE,
92    DI_SAVE_TABBED
93};
94
95enum LoadWhat { DI_LOAD_ALL, DI_LOAD_MARKED, DI_LOAD_LIST };
96
97class MatrixOrder : virtual Noncopyable {
98    GB_HASH *name2pos; // key = species name, value = order in sort_tree [1..n]
99                       // if no sort tree was specified, name2pos is NULL
100    int      leafs;    // number of leafs
101
102    bool tree_contains_dups; // unused (if true, matrix sorting works partly wrong)
103
104    void insert_in_hash(GBT_TREE *tree) {
105        if (tree->is_leaf) {
106            arb_assert(tree->name);
107            if (GBS_write_hash(name2pos, tree->name, ++leafs) != 0) {
108                tree_contains_dups = true;
109            }
110        }
111        else {
112            insert_in_hash(tree->rightson);
113            insert_in_hash(tree->leftson);
114        }
115    }
116
117public:
118    MatrixOrder(GBDATA *gb_main, GB_CSTR sort_tree_name);
119    ~MatrixOrder() { if (name2pos) GBS_free_hash(name2pos); }
120
121    bool defined() const { return leafs; }
122    int get_index(const char *name) const {
123        // return 1 for lowest and 'leafs' for highest species in sort-tee
124        // return 0 for all species missing in sort-tree
125        return defined() ? GBS_read_hash(name2pos, name) : -1;
126    }
127    void applyTo(struct TreeOrderedSpecies **gb_species_array, size_t array_size) const;
128};
129
130typedef void (*DI_MATRIX_CB)();
131
132class DI_MATRIX : virtual Noncopyable {
133    GBDATA  *gb_species_data;
134    long     seq_len;
135    char     cancel_columns[256];
136    size_t   entries_mem_size;
137    AliView *aliview;
138
139    GBDATA *get_gb_main() const { return aliview->get_gb_main(); }
140    double  corr(double dist, double b, double & sigma);
141    char   *calculate_overall_freqs(double rel_frequencies[AP_MAX], char *cancel_columns);
142    int     search_group(GBT_TREE *node, GB_HASH *hash, size_t& groupcnt, char *groupname, DI_ENTRY **groups);
143
144public:
145    // @@@ make members private:
146    bool             is_AA;
147    DI_ENTRY       **entries;
148    size_t           nentries;
149    AP_smatrix      *matrix;
150    DI_MATRIX_TYPE   matrix_type;
151
152    explicit DI_MATRIX(const AliView& aliview);
153    ~DI_MATRIX();
154
155    const char *get_aliname() const { return aliview->get_aliname(); }
156    const AliView *get_aliview() const { return aliview; }
157
158    GB_ERROR load(LoadWhat what, const MatrixOrder& order, bool show_warnings, GBDATA **species_list) __ATTR__USERESULT;
159    char *unload();
160    const char *save(const char *filename, enum DI_SAVE_TYPE type);
161
162    GB_ERROR  calculate(const char *cancel, DI_TRANSFORMATION transformation, bool *aborted_flag, AP_matrix *userdef_matrix);
163    GB_ERROR  calculate_pro(DI_TRANSFORMATION transformation, bool *aborted_flag);
164    GB_ERROR  extract_from_tree(const char *treename, bool *aborted_flag);
165
166    DI_TRANSFORMATION detect_transformation(std::string& msg);
167
168    char *compress(GBT_TREE *tree);
169};
170
171class DI_GLOBAL_MATRIX : virtual Noncopyable {
172    DI_MATRIX    *matrix;
173    DI_MATRIX_CB  changed_cb;
174
175    void announce_change() { if (changed_cb) changed_cb(); }
176
177    void forget_no_announce() {
178        delete matrix;
179        matrix = NULL;
180    }
181
182    void set(DI_MATRIX *new_global) { di_assert(!matrix); matrix = new_global; announce_change(); }
183
184public:
185    DI_GLOBAL_MATRIX() : matrix(NULL), changed_cb(NULL) {}
186    ~DI_GLOBAL_MATRIX() { forget(); }
187
188    DI_MATRIX *get() { return matrix; }
189    void forget() {
190        if (matrix) {
191            forget_no_announce();
192            announce_change();
193        }
194    }
195    void replaceBy(DI_MATRIX *new_global) { forget_no_announce(); set(new_global); }
196
197    bool exists() const { return matrix != NULL; }
198
199    void set_changed_cb(DI_MATRIX_CB cb) {
200        // announce_change(); // do by caller if really needed
201        changed_cb = cb;
202    }
203
204    DI_MATRIX *swap(DI_MATRIX *other) {
205        DI_MATRIX *prev = matrix;
206        matrix          = other;
207        announce_change();
208        return prev;
209    }
210
211    bool has_type(DI_MATRIX_TYPE type) const {
212        return matrix && matrix->matrix_type == type;
213    }
214    void forget_if_not_has_type(DI_MATRIX_TYPE wanted_type) {
215        if (matrix && matrix->matrix_type != wanted_type) {
216            forget();
217        }
218    }
219};
220
221extern DI_GLOBAL_MATRIX GLOBAL_MATRIX;
222
223class WeightedFilter;
224struct save_matrix_params {
225    const char           *awar_base;
226    const WeightedFilter *weighted_filter;
227};
228
229AW_window *DI_create_save_matrix_window(AW_root *aw_root, save_matrix_params *save_params);
230
231#else
232#error di_matr.hxx included twice
233#endif // DI_MATR_HXX
Note: See TracBrowser for help on using the repository browser.