| armadillo | | armadillo | |
| | | | |
| skipping to change at line 230 | | skipping to change at line 230 | |
| #include "armadillo_bits/op_cumsum_bones.hpp" | | #include "armadillo_bits/op_cumsum_bones.hpp" | |
| #include "armadillo_bits/op_symmat_bones.hpp" | | #include "armadillo_bits/op_symmat_bones.hpp" | |
| #include "armadillo_bits/op_hist_bones.hpp" | | #include "armadillo_bits/op_hist_bones.hpp" | |
| #include "armadillo_bits/op_unique_bones.hpp" | | #include "armadillo_bits/op_unique_bones.hpp" | |
| #include "armadillo_bits/op_toeplitz_bones.hpp" | | #include "armadillo_bits/op_toeplitz_bones.hpp" | |
| #include "armadillo_bits/op_fft_bones.hpp" | | #include "armadillo_bits/op_fft_bones.hpp" | |
| #include "armadillo_bits/op_any_bones.hpp" | | #include "armadillo_bits/op_any_bones.hpp" | |
| #include "armadillo_bits/op_all_bones.hpp" | | #include "armadillo_bits/op_all_bones.hpp" | |
| #include "armadillo_bits/op_normalise_bones.hpp" | | #include "armadillo_bits/op_normalise_bones.hpp" | |
| #include "armadillo_bits/op_clamp_bones.hpp" | | #include "armadillo_bits/op_clamp_bones.hpp" | |
|
| | | #include "armadillo_bits/op_expmat_bones.hpp" | |
| | | | |
| #include "armadillo_bits/glue_times_bones.hpp" | | #include "armadillo_bits/glue_times_bones.hpp" | |
| #include "armadillo_bits/glue_mixed_bones.hpp" | | #include "armadillo_bits/glue_mixed_bones.hpp" | |
| #include "armadillo_bits/glue_cov_bones.hpp" | | #include "armadillo_bits/glue_cov_bones.hpp" | |
| #include "armadillo_bits/glue_cor_bones.hpp" | | #include "armadillo_bits/glue_cor_bones.hpp" | |
| #include "armadillo_bits/glue_kron_bones.hpp" | | #include "armadillo_bits/glue_kron_bones.hpp" | |
| #include "armadillo_bits/glue_cross_bones.hpp" | | #include "armadillo_bits/glue_cross_bones.hpp" | |
| #include "armadillo_bits/glue_join_bones.hpp" | | #include "armadillo_bits/glue_join_bones.hpp" | |
| #include "armadillo_bits/glue_relational_bones.hpp" | | #include "armadillo_bits/glue_relational_bones.hpp" | |
| #include "armadillo_bits/glue_solve_bones.hpp" | | #include "armadillo_bits/glue_solve_bones.hpp" | |
| | | | |
| skipping to change at line 430 | | skipping to change at line 431 | |
| #include "armadillo_bits/fn_any.hpp" | | #include "armadillo_bits/fn_any.hpp" | |
| #include "armadillo_bits/fn_all.hpp" | | #include "armadillo_bits/fn_all.hpp" | |
| #include "armadillo_bits/fn_size.hpp" | | #include "armadillo_bits/fn_size.hpp" | |
| #include "armadillo_bits/fn_numel.hpp" | | #include "armadillo_bits/fn_numel.hpp" | |
| #include "armadillo_bits/fn_inplace_strans.hpp" | | #include "armadillo_bits/fn_inplace_strans.hpp" | |
| #include "armadillo_bits/fn_inplace_trans.hpp" | | #include "armadillo_bits/fn_inplace_trans.hpp" | |
| #include "armadillo_bits/fn_randi.hpp" | | #include "armadillo_bits/fn_randi.hpp" | |
| #include "armadillo_bits/fn_cond.hpp" | | #include "armadillo_bits/fn_cond.hpp" | |
| #include "armadillo_bits/fn_normalise.hpp" | | #include "armadillo_bits/fn_normalise.hpp" | |
| #include "armadillo_bits/fn_clamp.hpp" | | #include "armadillo_bits/fn_clamp.hpp" | |
|
| | | #include "armadillo_bits/fn_expmat.hpp" | |
| | | | |
| #include "armadillo_bits/fn_speye.hpp" | | #include "armadillo_bits/fn_speye.hpp" | |
| #include "armadillo_bits/fn_spones.hpp" | | #include "armadillo_bits/fn_spones.hpp" | |
| #include "armadillo_bits/fn_sprandn.hpp" | | #include "armadillo_bits/fn_sprandn.hpp" | |
| #include "armadillo_bits/fn_sprandu.hpp" | | #include "armadillo_bits/fn_sprandu.hpp" | |
| #include "armadillo_bits/fn_eigs_sym.hpp" | | #include "armadillo_bits/fn_eigs_sym.hpp" | |
| #include "armadillo_bits/fn_eigs_gen.hpp" | | #include "armadillo_bits/fn_eigs_gen.hpp" | |
| #include "armadillo_bits/fn_norm_sparse.hpp" | | #include "armadillo_bits/fn_norm_sparse.hpp" | |
| | | | |
| // | | // | |
| | | | |
| skipping to change at line 545 | | skipping to change at line 547 | |
| #include "armadillo_bits/op_cumsum_meat.hpp" | | #include "armadillo_bits/op_cumsum_meat.hpp" | |
| #include "armadillo_bits/op_symmat_meat.hpp" | | #include "armadillo_bits/op_symmat_meat.hpp" | |
| #include "armadillo_bits/op_hist_meat.hpp" | | #include "armadillo_bits/op_hist_meat.hpp" | |
| #include "armadillo_bits/op_unique_meat.hpp" | | #include "armadillo_bits/op_unique_meat.hpp" | |
| #include "armadillo_bits/op_toeplitz_meat.hpp" | | #include "armadillo_bits/op_toeplitz_meat.hpp" | |
| #include "armadillo_bits/op_fft_meat.hpp" | | #include "armadillo_bits/op_fft_meat.hpp" | |
| #include "armadillo_bits/op_any_meat.hpp" | | #include "armadillo_bits/op_any_meat.hpp" | |
| #include "armadillo_bits/op_all_meat.hpp" | | #include "armadillo_bits/op_all_meat.hpp" | |
| #include "armadillo_bits/op_normalise_meat.hpp" | | #include "armadillo_bits/op_normalise_meat.hpp" | |
| #include "armadillo_bits/op_clamp_meat.hpp" | | #include "armadillo_bits/op_clamp_meat.hpp" | |
|
| | | #include "armadillo_bits/op_expmat_meat.hpp" | |
| | | | |
| #include "armadillo_bits/glue_times_meat.hpp" | | #include "armadillo_bits/glue_times_meat.hpp" | |
| #include "armadillo_bits/glue_mixed_meat.hpp" | | #include "armadillo_bits/glue_mixed_meat.hpp" | |
| #include "armadillo_bits/glue_cov_meat.hpp" | | #include "armadillo_bits/glue_cov_meat.hpp" | |
| #include "armadillo_bits/glue_cor_meat.hpp" | | #include "armadillo_bits/glue_cor_meat.hpp" | |
| #include "armadillo_bits/glue_kron_meat.hpp" | | #include "armadillo_bits/glue_kron_meat.hpp" | |
| #include "armadillo_bits/glue_cross_meat.hpp" | | #include "armadillo_bits/glue_cross_meat.hpp" | |
| #include "armadillo_bits/glue_join_meat.hpp" | | #include "armadillo_bits/glue_join_meat.hpp" | |
| #include "armadillo_bits/glue_relational_meat.hpp" | | #include "armadillo_bits/glue_relational_meat.hpp" | |
| #include "armadillo_bits/glue_solve_meat.hpp" | | #include "armadillo_bits/glue_solve_meat.hpp" | |
| | | | |
End of changes. 3 change blocks. |
| 0 lines changed or deleted | | 3 lines changed or added | |
|
| gmm_diag_bones.hpp | | gmm_diag_bones.hpp | |
| | | | |
| skipping to change at line 46 | | skipping to change at line 46 | |
| static const gmm_seed_random_spread random_spread; | | static const gmm_seed_random_spread random_spread; | |
| | | | |
| namespace gmm_priv | | namespace gmm_priv | |
| { | | { | |
| | | | |
| struct gmm_empty_arg {}; | | struct gmm_empty_arg {}; | |
| | | | |
| #if defined(_OPENMP) | | #if defined(_OPENMP) | |
| struct arma_omp_state | | struct arma_omp_state | |
| { | | { | |
|
| const int dynamic_state; | | const int orig_dynamic_state; | |
| | | | |
|
| inline arma_omp_state() : dynamic_state(omp_get_dynamic()) {} | | inline arma_omp_state() : orig_dynamic_state(omp_get_dynamic()) { omp_ | |
| inline ~arma_omp_state() { omp_set_dynamic(dynamic_state); } | | set_dynamic(0); } | |
| | | inline ~arma_omp_state() { omp_set_dynamic(orig_dynamic_state); } | |
| }; | | }; | |
|
| | | #else | |
| | | struct arma_omp_state {}; | |
| #endif | | #endif | |
| | | | |
| template<typename eT> | | template<typename eT> | |
| class gmm_diag | | class gmm_diag | |
| { | | { | |
| public: | | public: | |
| | | | |
| arma_aligned const Mat<eT> means; | | arma_aligned const Mat<eT> means; | |
| arma_aligned const Mat<eT> dcovs; | | arma_aligned const Mat<eT> dcovs; | |
| arma_aligned const Row<eT> hefts; | | arma_aligned const Row<eT> hefts; | |
| | | | |
| skipping to change at line 138 | | skipping to change at line 139 | |
| arma_aligned Col<eT> mah_aux; | | arma_aligned Col<eT> mah_aux; | |
| | | | |
| // | | // | |
| | | | |
| inline void init(const gmm_diag& x); | | inline void init(const gmm_diag& x); | |
| | | | |
| inline void init(const uword in_n_dim, const uword in_n_gaus); | | inline void init(const uword in_n_dim, const uword in_n_gaus); | |
| | | | |
| inline void init_constants(); | | inline void init_constants(); | |
| | | | |
|
| | | inline umat internal_gen_boundaries(const uword N) const; | |
| | | | |
| inline eT internal_scalar_log_p(const eT* x ) const; | | inline eT internal_scalar_log_p(const eT* x ) const; | |
| inline eT internal_scalar_log_p(const eT* x, const uword gaus_id) const; | | inline eT internal_scalar_log_p(const eT* x, const uword gaus_id) const; | |
| | | | |
| template<typename T1> inline Row<eT> internal_vec_log_p(const T1& X
) const; | | template<typename T1> inline Row<eT> internal_vec_log_p(const T1& X
) const; | |
| template<typename T1> inline Row<eT> internal_vec_log_p(const T1& X, cons
t uword gaus_id) const; | | template<typename T1> inline Row<eT> internal_vec_log_p(const T1& X, cons
t uword gaus_id) const; | |
| | | | |
| template<typename T1> inline eT internal_avg_log_p(const T1& X
) const; | | template<typename T1> inline eT internal_avg_log_p(const T1& X
) const; | |
| template<typename T1> inline eT internal_avg_log_p(const T1& X, const uwo
rd gaus_id) const; | | template<typename T1> inline eT internal_avg_log_p(const T1& X, const uwo
rd gaus_id) const; | |
| | | | |
| template<typename T1> inline uword internal_scalar_assign(const T1& X, co
nst gmm_dist_mode& dist_mode) const; | | template<typename T1> inline uword internal_scalar_assign(const T1& X, co
nst gmm_dist_mode& dist_mode) const; | |
| | | | |
| skipping to change at line 167 | | skipping to change at line 170 | |
| template<uword dist_id> inline void generate_initial_dcovs_and_hefts(cons
t Mat<eT>& X, const eT var_floor); | | template<uword dist_id> inline void generate_initial_dcovs_and_hefts(cons
t Mat<eT>& X, const eT var_floor); | |
| | | | |
| template<uword dist_id> inline bool km_iterate(const Mat<eT>& X, const uw
ord max_iter, const bool verbose); | | template<uword dist_id> inline bool km_iterate(const Mat<eT>& X, const uw
ord max_iter, const bool verbose); | |
| | | | |
| template<uword dist_id> inline void km_update_stats(const Mat<eT>& X, con
st uword start_index, const uword end_index, const Mat<eT>& old_means, fiel
d< running_mean_vec<eT> >& running_means) const; | | template<uword dist_id> inline void km_update_stats(const Mat<eT>& X, con
st uword start_index, const uword end_index, const Mat<eT>& old_means, fiel
d< running_mean_vec<eT> >& running_means) const; | |
| | | | |
| // | | // | |
| | | | |
| inline bool em_iterate(const Mat<eT>& X, const uword max_iter, const eT v
ar_floor, const bool verbose); | | inline bool em_iterate(const Mat<eT>& X, const uword max_iter, const eT v
ar_floor, const bool verbose); | |
| | | | |
|
| inline void em_update_params(const Mat<eT>& X, const field<uvec>& t_bound
ary, field< Mat<eT> >& t_acc_means, field< Mat<eT> >& t_acc_dcovs, field< C
ol<eT> >& t_acc_norm_lhoods, field< Col<eT> >& t_gaus_log_lhoods, Col<eT>&
t_progress_log_lhoods); | | inline void em_update_params(const Mat<eT>& X, const umat& boundaries, fi
eld< Mat<eT> >& t_acc_means, field< Mat<eT> >& t_acc_dcovs, field< Col<eT>
>& t_acc_norm_lhoods, field< Col<eT> >& t_gaus_log_lhoods, Col<eT>& t_progr
ess_log_lhoods); | |
| | | | |
|
| inline void em_generate_acc(const Mat<eT>& X, const uvec& boundary, Mat<e
T>& acc_means, Mat<eT>& acc_dcovs, Col<eT>& acc_norm_lhoods, Col<eT>& gaus_
log_lhoods, eT& progress_log_lhood) const; | | inline void em_generate_acc(const Mat<eT>& X, const uword start_index, co
nst uword end_index, Mat<eT>& acc_means, Mat<eT>& acc_dcovs, Col<eT>& acc_n
orm_lhoods, Col<eT>& gaus_log_lhoods, eT& progress_log_lhood) const; | |
| | | | |
| inline void em_fix_params(const eT var_floor); | | inline void em_fix_params(const eT var_floor); | |
| }; | | }; | |
| | | | |
| } | | } | |
| | | | |
| typedef gmm_priv::gmm_diag<double> gmm_diag; | | typedef gmm_priv::gmm_diag<double> gmm_diag; | |
| typedef gmm_priv::gmm_diag<float> fgmm_diag; | | typedef gmm_priv::gmm_diag<float> fgmm_diag; | |
| | | | |
End of changes. 6 change blocks. |
| 6 lines changed or deleted | | 10 lines changed or added | |
|
| gmm_diag_meat.hpp | | gmm_diag_meat.hpp | |
| | | | |
| skipping to change at line 759 | | skipping to change at line 759 | |
| { | | { | |
| const eT logdet = accu( log(dcovs.col(i)) ); | | const eT logdet = accu( log(dcovs.col(i)) ); | |
| | | | |
| log_det_etc[i] = eT(-1) * ( tmp + eT(0.5) * logdet ); | | log_det_etc[i] = eT(-1) * ( tmp + eT(0.5) * logdet ); | |
| } | | } | |
| | | | |
| log_hefts = log(hefts); // TODO: possible issue when one of the hefts is
zero | | log_hefts = log(hefts); // TODO: possible issue when one of the hefts is
zero | |
| } | | } | |
| | | | |
| template<typename eT> | | template<typename eT> | |
|
| | | inline | |
| | | umat | |
| | | gmm_diag<eT>::internal_gen_boundaries(const uword N) const | |
| | | { | |
| | | arma_extra_debug_sigprint(); | |
| | | | |
| | | #if defined(_OPENMP) | |
| | | // const uword n_cores = 0; | |
| | | const uword n_cores = uword(omp_get_num_procs()); | |
| | | const uword n_threads = (n_cores > 0) ? ( (n_cores <= N) ? n_cores : 1 | |
| | | ) : 1; | |
| | | #else | |
| | | // static const uword n_cores = 0; | |
| | | static const uword n_threads = 1; | |
| | | #endif | |
| | | | |
| | | // get_stream_err2() << "gmm_diag::internal_gen_boundaries(): n_cores: | |
| | | " << n_cores << '\n'; | |
| | | // get_stream_err2() << "gmm_diag::internal_gen_boundaries(): n_threads: | |
| | | " << n_threads << '\n'; | |
| | | | |
| | | umat boundaries(2, n_threads); | |
| | | | |
| | | if(N > 0) | |
| | | { | |
| | | const uword chunk_size = N / n_threads; | |
| | | | |
| | | uword count = 0; | |
| | | | |
| | | for(uword t=0; t<n_threads; t++) | |
| | | { | |
| | | boundaries.at(0,t) = count; | |
| | | | |
| | | count += chunk_size; | |
| | | | |
| | | boundaries.at(1,t) = count-1; | |
| | | } | |
| | | | |
| | | boundaries.at(1,n_threads-1) = N - 1; | |
| | | } | |
| | | else | |
| | | { | |
| | | boundaries.zeros(); | |
| | | } | |
| | | | |
| | | // get_stream_err2() << "gmm_diag::internal_gen_boundaries(): boundaries: | |
| | | " << '\n' << boundaries << '\n'; | |
| | | | |
| | | return boundaries; | |
| | | } | |
| | | | |
| | | template<typename eT> | |
| arma_hot | | arma_hot | |
| inline | | inline | |
| eT | | eT | |
| gmm_diag<eT>::internal_scalar_log_p(const eT* x) const | | gmm_diag<eT>::internal_scalar_log_p(const eT* x) const | |
| { | | { | |
| arma_extra_debug_sigprint(); | | arma_extra_debug_sigprint(); | |
| | | | |
| const eT* log_hefts_mem = log_hefts.mem; | | const eT* log_hefts_mem = log_hefts.mem; | |
| | | | |
| const uword N_gaus = means.n_cols; | | const uword N_gaus = means.n_cols; | |
| | | | |
| skipping to change at line 843 | | skipping to change at line 891 | |
| gmm_diag<eT>::internal_vec_log_p(const T1& X) const | | gmm_diag<eT>::internal_vec_log_p(const T1& X) const | |
| { | | { | |
| arma_extra_debug_sigprint(); | | arma_extra_debug_sigprint(); | |
| | | | |
| arma_debug_check( (X.n_rows != means.n_rows), "gmm_diag::log_p(): incompa
tible dimensions" ); | | arma_debug_check( (X.n_rows != means.n_rows), "gmm_diag::log_p(): incompa
tible dimensions" ); | |
| | | | |
| const uword N = X.n_cols; | | const uword N = X.n_cols; | |
| | | | |
| Row<eT> out(N); | | Row<eT> out(N); | |
| | | | |
|
| eT* out_mem = out.memptr(); | | if(N > 0) | |
| | | | |
| for(uword i=0; i < N; ++i) | | | |
| { | | { | |
|
| out_mem[i] = internal_scalar_log_p( X.colptr(i) ); | | #if defined(_OPENMP) | |
| | | { | |
| | | const arma_omp_state save_omp_state; | |
| | | | |
| | | const umat boundaries = internal_gen_boundaries(N); | |
| | | | |
| | | const uword n_threads = boundaries.n_cols; | |
| | | | |
| | | #pragma omp parallel for | |
| | | for(uword t=0; t < n_threads; ++t) | |
| | | { | |
| | | const uword start_index = boundaries.at(0,t); | |
| | | const uword end_index = boundaries.at(1,t); | |
| | | | |
| | | eT* out_mem = out.memptr(); | |
| | | | |
| | | for(uword i=start_index; i <= end_index; ++i) | |
| | | { | |
| | | out_mem[i] = internal_scalar_log_p( X.colptr(i) ); | |
| | | } | |
| | | } | |
| | | } | |
| | | #else | |
| | | { | |
| | | eT* out_mem = out.memptr(); | |
| | | | |
| | | for(uword i=0; i < N; ++i) | |
| | | { | |
| | | out_mem[i] = internal_scalar_log_p( X.colptr(i) ); | |
| | | } | |
| | | } | |
| | | #endif | |
| } | | } | |
| | | | |
| return out; | | return out; | |
| } | | } | |
| | | | |
| template<typename eT> | | template<typename eT> | |
| template<typename T1> | | template<typename T1> | |
| inline | | inline | |
| Row<eT> | | Row<eT> | |
| gmm_diag<eT>::internal_vec_log_p(const T1& X, const uword gaus_id) const | | gmm_diag<eT>::internal_vec_log_p(const T1& X, const uword gaus_id) const | |
| { | | { | |
| arma_extra_debug_sigprint(); | | arma_extra_debug_sigprint(); | |
| | | | |
| arma_debug_check( (X.n_rows != means.n_rows), "gmm_diag::log_p(): incompa
tible dimensions" ); | | arma_debug_check( (X.n_rows != means.n_rows), "gmm_diag::log_p(): incompa
tible dimensions" ); | |
| arma_debug_check( (gaus_id >= means.n_cols), "gmm_diag::log_p(): gaus_id
is out of range" ); | | arma_debug_check( (gaus_id >= means.n_cols), "gmm_diag::log_p(): gaus_id
is out of range" ); | |
| | | | |
| const uword N = X.n_cols; | | const uword N = X.n_cols; | |
| | | | |
| Row<eT> out(N); | | Row<eT> out(N); | |
| | | | |
|
| eT* out_mem = out.memptr(); | | if(N > 0) | |
| | | | |
| for(uword i=0; i < N; ++i) | | | |
| { | | { | |
|
| out_mem[i] = internal_scalar_log_p( X.colptr(i), gaus_id ); | | #if defined(_OPENMP) | |
| | | { | |
| | | const arma_omp_state save_omp_state; | |
| | | | |
| | | const umat boundaries = internal_gen_boundaries(N); | |
| | | | |
| | | const uword n_threads = boundaries.n_cols; | |
| | | | |
| | | #pragma omp parallel for | |
| | | for(uword t=0; t < n_threads; ++t) | |
| | | { | |
| | | const uword start_index = boundaries.at(0,t); | |
| | | const uword end_index = boundaries.at(1,t); | |
| | | | |
| | | eT* out_mem = out.memptr(); | |
| | | | |
| | | for(uword i=start_index; i <= end_index; ++i) | |
| | | { | |
| | | out_mem[i] = internal_scalar_log_p( X.colptr(i), gaus_id ); | |
| | | } | |
| | | } | |
| | | } | |
| | | #else | |
| | | { | |
| | | eT* out_mem = out.memptr(); | |
| | | | |
| | | for(uword i=0; i < N; ++i) | |
| | | { | |
| | | out_mem[i] = internal_scalar_log_p( X.colptr(i), gaus_id ); | |
| | | } | |
| | | } | |
| | | #endif | |
| } | | } | |
| | | | |
| return out; | | return out; | |
| } | | } | |
| | | | |
| template<typename eT> | | template<typename eT> | |
| template<typename T1> | | template<typename T1> | |
| inline | | inline | |
| eT | | eT | |
| gmm_diag<eT>::internal_avg_log_p(const T1& X) const | | gmm_diag<eT>::internal_avg_log_p(const T1& X) const | |
| { | | { | |
| arma_extra_debug_sigprint(); | | arma_extra_debug_sigprint(); | |
| | | | |
| arma_debug_check( (X.n_rows != means.n_rows), "gmm_diag::avg_log_p(): inc
ompatible dimensions" ); | | arma_debug_check( (X.n_rows != means.n_rows), "gmm_diag::avg_log_p(): inc
ompatible dimensions" ); | |
| | | | |
|
| running_mean_scalar<eT> running_mean; | | | |
| | | | |
| const uword N = X.n_cols; | | const uword N = X.n_cols; | |
| | | | |
|
| for(uword i=0; i<N; ++i) | | if(N == 0) { return (-Datum<eT>::inf); } | |
| | | | |
| | | #if defined(_OPENMP) | |
| { | | { | |
|
| running_mean( internal_scalar_log_p( X.colptr(i) ) ); | | const arma_omp_state save_omp_state; | |
| | | | |
| | | const umat boundaries = internal_gen_boundaries(N); | |
| | | | |
| | | const uword n_threads = boundaries.n_cols; | |
| | | | |
| | | field< running_mean_scalar<eT> > t_running_means(n_threads); | |
| | | | |
| | | #pragma omp parallel for | |
| | | for(uword t=0; t < n_threads; ++t) | |
| | | { | |
| | | const uword start_index = boundaries.at(0,t); | |
| | | const uword end_index = boundaries.at(1,t); | |
| | | | |
| | | running_mean_scalar<eT>& current_running_mean = t_running_means[t]; | |
| | | | |
| | | for(uword i=start_index; i <= end_index; ++i) | |
| | | { | |
| | | current_running_mean( internal_scalar_log_p( X.colptr(i) ) ); | |
| | | } | |
| | | } | |
| | | | |
| | | eT avg = eT(0); | |
| | | | |
| | | for(uword t=0; t < n_threads; ++t) | |
| | | { | |
| | | running_mean_scalar<eT>& current_running_mean = t_running_means[t]; | |
| | | | |
| | | const eT w = eT(current_running_mean.count()) / eT(N); | |
| | | | |
| | | avg += w * current_running_mean.mean(); | |
| | | } | |
| | | | |
| | | return avg; | |
| } | | } | |
|
| | | #else | |
| | | { | |
| | | running_mean_scalar<eT> running_mean; | |
| | | | |
|
| return ( (N > 0) ? running_mean.mean() : (-Datum<eT>::inf) ); | | for(uword i=0; i<N; ++i) | |
| | | { | |
| | | running_mean( internal_scalar_log_p( X.colptr(i) ) ); | |
| | | } | |
| | | | |
| | | return running_mean.mean(); | |
| | | } | |
| | | #endif | |
| } | | } | |
| | | | |
| template<typename eT> | | template<typename eT> | |
| template<typename T1> | | template<typename T1> | |
| inline | | inline | |
| eT | | eT | |
| gmm_diag<eT>::internal_avg_log_p(const T1& X, const uword gaus_id) const | | gmm_diag<eT>::internal_avg_log_p(const T1& X, const uword gaus_id) const | |
| { | | { | |
| arma_extra_debug_sigprint(); | | arma_extra_debug_sigprint(); | |
| | | | |
| arma_debug_check( (X.n_rows != means.n_rows), "gmm_diag::avg_log_p(): inc
ompatible dimensions" ); | | arma_debug_check( (X.n_rows != means.n_rows), "gmm_diag::avg_log_p(): inc
ompatible dimensions" ); | |
| arma_debug_check( (gaus_id >= means.n_cols), "gmm_diag::avg_log_p(): spe
cified gaussian is out of range" ); | | arma_debug_check( (gaus_id >= means.n_cols), "gmm_diag::avg_log_p(): spe
cified gaussian is out of range" ); | |
| | | | |
|
| running_mean_scalar<eT> running_mean; | | | |
| | | | |
| const uword N = X.n_cols; | | const uword N = X.n_cols; | |
| | | | |
|
| for(uword i=0; i<N; ++i) | | if(N == 0) { return (-Datum<eT>::inf); } | |
| | | | |
| | | #if defined(_OPENMP) | |
| { | | { | |
|
| running_mean( internal_scalar_log_p( X.colptr(i), gaus_id ) ); | | const arma_omp_state save_omp_state; | |
| | | | |
| | | const umat boundaries = internal_gen_boundaries(N); | |
| | | | |
| | | const uword n_threads = boundaries.n_cols; | |
| | | | |
| | | field< running_mean_scalar<eT> > t_running_means(n_threads); | |
| | | | |
| | | #pragma omp parallel for | |
| | | for(uword t=0; t < n_threads; ++t) | |
| | | { | |
| | | const uword start_index = boundaries.at(0,t); | |
| | | const uword end_index = boundaries.at(1,t); | |
| | | | |
| | | running_mean_scalar<eT>& current_running_mean = t_running_means[t]; | |
| | | | |
| | | for(uword i=start_index; i <= end_index; ++i) | |
| | | { | |
| | | current_running_mean( internal_scalar_log_p( X.colptr(i), gaus_id) | |
| | | ); | |
| | | } | |
| | | } | |
| | | | |
| | | eT avg = eT(0); | |
| | | | |
| | | for(uword t=0; t < n_threads; ++t) | |
| | | { | |
| | | running_mean_scalar<eT>& current_running_mean = t_running_means[t]; | |
| | | | |
| | | const eT w = eT(current_running_mean.count()) / eT(N); | |
| | | | |
| | | avg += w * current_running_mean.mean(); | |
| | | } | |
| | | | |
| | | return avg; | |
| } | | } | |
|
| | | #else | |
| | | { | |
| | | running_mean_scalar<eT> running_mean; | |
| | | | |
|
| return ( (N > 0) ? running_mean.mean() : (-Datum<eT>::inf) ); | | const uword N = X.n_cols; | |
| | | | |
| | | for(uword i=0; i<N; ++i) | |
| | | { | |
| | | running_mean( internal_scalar_log_p( X.colptr(i), gaus_id ) ); | |
| | | } | |
| | | | |
| | | return running_mean.mean(); | |
| | | } | |
| | | #endif | |
| } | | } | |
| | | | |
| template<typename eT> | | template<typename eT> | |
| template<typename T1> | | template<typename T1> | |
| inline | | inline | |
| uword | | uword | |
| gmm_diag<eT>::internal_scalar_assign(const T1& X, const gmm_dist_mode& dist
_mode) const | | gmm_diag<eT>::internal_scalar_assign(const T1& X, const gmm_dist_mode& dist
_mode) const | |
| { | | { | |
| arma_extra_debug_sigprint(); | | arma_extra_debug_sigprint(); | |
| | | | |
| | | | |
| skipping to change at line 1285 | | skipping to change at line 1479 | |
| | | | |
| running_mean_scalar<double> rs_delta; | | running_mean_scalar<double> rs_delta; | |
| | | | |
| field< running_mean_vec<eT> > running_means(N_gaus); | | field< running_mean_vec<eT> > running_means(N_gaus); | |
| | | | |
| const eT* mah_aux_mem = mah_aux.memptr(); | | const eT* mah_aux_mem = mah_aux.memptr(); | |
| | | | |
| #if defined(_OPENMP) | | #if defined(_OPENMP) | |
| const arma_omp_state save_omp_state; | | const arma_omp_state save_omp_state; | |
| | | | |
|
| omp_set_dynamic(0); | | const umat boundaries = internal_gen_boundaries(X.n_cols); | |
| | | | |
|
| //const uword n_cores = 0; | | const uword n_threads = boundaries.n_cols; | |
| const uword n_cores = uword(omp_get_num_procs()); | | | |
| const uword n_threads = (n_cores > 0) ? ( (n_cores <= X.n_cols) ? n_cor | | | |
| es : 1 ) : 1; | | | |
| | | | |
| field< field< running_mean_vec<eT> > > t_running_means(n_threads); | | field< field< running_mean_vec<eT> > > t_running_means(n_threads); | |
| | | | |
| for(uword t=0; t < n_threads; ++t) { t_running_means[t].set_size(N_gau
s); } | | for(uword t=0; t < n_threads; ++t) { t_running_means[t].set_size(N_gau
s); } | |
| | | | |
|
| field< uvec > t_boundary(n_threads); | | | |
| | | | |
| const uword chunk_size = X.n_cols / n_threads; | | | |
| | | | |
| uword vec_count = 0; | | | |
| | | | |
| for(uword t=0; t<n_threads; t++) | | | |
| { | | | |
| t_boundary[t].set_size(2); | | | |
| | | | |
| t_boundary[t][0] = vec_count; | | | |
| | | | |
| vec_count += chunk_size; | | | |
| | | | |
| t_boundary[t][1] = vec_count-1; | | | |
| } | | | |
| | | | |
| t_boundary[n_threads-1][1] = X.n_cols - 1; | | | |
| | | | |
| vec tmp_mean(N_dims); | | vec tmp_mean(N_dims); | |
| | | | |
| if(verbose) | | if(verbose) | |
| { | | { | |
|
| get_stream_err2() << "gmm_diag::learn(): k-means: n_threads: " << n_ | | get_stream_err2() << "gmm_diag::learn(): k-means: n_threads: " << n_t | |
| threads << '\n'; | | hreads << '\n'; | |
| get_stream_err2() << "gmm_diag::learn(): k-means: chunk_size: " << ch | | | |
| unk_size << '\n'; | | | |
| } | | } | |
| #endif | | #endif | |
| | | | |
| for(uword iter=1; iter <= max_iter; ++iter) | | for(uword iter=1; iter <= max_iter; ++iter) | |
| { | | { | |
| #if defined(_OPENMP) | | #if defined(_OPENMP) | |
| { | | { | |
| for(uword t=0; t < n_threads; ++t) | | for(uword t=0; t < n_threads; ++t) | |
| { | | { | |
| for(uword g=0; g < N_gaus; ++g) { t_running_means[t][g].reset(); } | | for(uword g=0; g < N_gaus; ++g) { t_running_means[t][g].reset(); } | |
| } | | } | |
| | | | |
| // km_update_stats() is the "map" operation, which produces partial m
eans | | // km_update_stats() is the "map" operation, which produces partial m
eans | |
| | | | |
| #pragma omp parallel for | | #pragma omp parallel for | |
| for(uword t=0; t < n_threads; ++t) | | for(uword t=0; t < n_threads; ++t) | |
| { | | { | |
|
| const uvec& boundary = t_boundary[t]; | | | |
| | | | |
| field< running_mean_vec<eT> >& current_running_means = t_running_me
ans[t]; | | field< running_mean_vec<eT> >& current_running_means = t_running_me
ans[t]; | |
| | | | |
|
| km_update_stats<dist_id>(X, boundary[0], boundary[1], old_means, cu
rrent_running_means); | | km_update_stats<dist_id>(X, boundaries.at(0,t), boundaries.at(1,t),
old_means, current_running_means); | |
| } | | } | |
| | | | |
| // the "reduce" operation, which combines the partial means produced
by the separate threads; | | // the "reduce" operation, which combines the partial means produced
by the separate threads; | |
| // takes into account the counts for each mean | | // takes into account the counts for each mean | |
| | | | |
| for(uword g=0; g < N_gaus; ++g) | | for(uword g=0; g < N_gaus; ++g) | |
| { | | { | |
| uword total_count = 0; | | uword total_count = 0; | |
| | | | |
| for(uword t=0; t < n_threads; ++t) { total_count += t_running_mean
s[t][g].count(); } | | for(uword t=0; t < n_threads; ++t) { total_count += t_running_mean
s[t][g].count(); } | |
| | | | |
| skipping to change at line 1544 | | skipping to change at line 1714 | |
| get_stream_err2().unsetf(ios::uppercase); | | get_stream_err2().unsetf(ios::uppercase); | |
| get_stream_err2().unsetf(ios::showpos); | | get_stream_err2().unsetf(ios::showpos); | |
| get_stream_err2().unsetf(ios::scientific); | | get_stream_err2().unsetf(ios::scientific); | |
| | | | |
| get_stream_err2().setf(ios::right); | | get_stream_err2().setf(ios::right); | |
| get_stream_err2().setf(ios::fixed); | | get_stream_err2().setf(ios::fixed); | |
| } | | } | |
| | | | |
| #if defined(_OPENMP) | | #if defined(_OPENMP) | |
| const arma_omp_state save_omp_state; | | const arma_omp_state save_omp_state; | |
|
| | | | |
| omp_set_dynamic(0); | | | |
| | | | |
| // const uword n_cores = 0; | | | |
| const uword n_cores = uword(omp_get_num_procs()); | | | |
| const uword n_threads = (n_cores > 0) ? ( (n_cores <= X.n_cols) ? n_cor | | | |
| es : 1 ) : 1; | | | |
| #else | | | |
| // static const uword n_cores = 0; | | | |
| static const uword n_threads = 1; | | | |
| #endif | | #endif | |
| | | | |
|
| // get_stream_err2() << "n_cores: " << n_cores << '\n'; | | const umat boundaries = internal_gen_boundaries(X.n_cols); | |
| // get_stream_err2() << "n_threads: " << n_threads << '\n'; | | | |
| | | | |
|
| field< uvec > t_boundary(n_threads); | | const uword n_threads = boundaries.n_cols; | |
| | | | |
| field< Mat<eT> > t_acc_means(n_threads); | | field< Mat<eT> > t_acc_means(n_threads); | |
| field< Mat<eT> > t_acc_dcovs(n_threads); | | field< Mat<eT> > t_acc_dcovs(n_threads); | |
| | | | |
| field< Col<eT> > t_acc_norm_lhoods(n_threads); | | field< Col<eT> > t_acc_norm_lhoods(n_threads); | |
| field< Col<eT> > t_gaus_log_lhoods(n_threads); | | field< Col<eT> > t_gaus_log_lhoods(n_threads); | |
| | | | |
| Col<eT> t_progress_log_lhood(n_threads); | | Col<eT> t_progress_log_lhood(n_threads); | |
| | | | |
| for(uword t=0; t<n_threads; t++) | | for(uword t=0; t<n_threads; t++) | |
| { | | { | |
|
| t_boundary[t].set_size(2); | | | |
| | | | |
| t_acc_means[t].set_size(N_dims, N_gaus); | | t_acc_means[t].set_size(N_dims, N_gaus); | |
| t_acc_dcovs[t].set_size(N_dims, N_gaus); | | t_acc_dcovs[t].set_size(N_dims, N_gaus); | |
| | | | |
| t_acc_norm_lhoods[t].set_size(N_gaus); | | t_acc_norm_lhoods[t].set_size(N_gaus); | |
| t_gaus_log_lhoods[t].set_size(N_gaus); | | t_gaus_log_lhoods[t].set_size(N_gaus); | |
| } | | } | |
| | | | |
|
| const uword chunk_size = X.n_cols / n_threads; | | | |
| | | | |
| uword count = 0; | | | |
| | | | |
| for(uword t=0; t<n_threads; t++) | | | |
| { | | | |
| t_boundary[t][0] = count; | | | |
| | | | |
| count += chunk_size; | | | |
| | | | |
| t_boundary[t][1] = count-1; | | | |
| } | | | |
| | | | |
| t_boundary[n_threads-1][1] = X.n_cols - 1; | | | |
| | | | |
| // get_stream_err2() << "t_boundary.n_elem: " << t_boundary.n_elem << '\n | | | |
| '; | | | |
| // t_boundary.print("t_boundary:"); | | | |
| | | | |
| #if defined(_OPENMP) | | #if defined(_OPENMP) | |
| if(verbose) | | if(verbose) | |
| { | | { | |
|
| get_stream_err2() << "gmm_diag::learn(): EM: n_threads: " << n_threa | | get_stream_err2() << "gmm_diag::learn(): EM: n_threads: " << n_thread | |
| ds << '\n'; | | s << '\n'; | |
| get_stream_err2() << "gmm_diag::learn(): EM: chunk_size: " << chunk_s | | | |
| ize << '\n'; | | | |
| } | | } | |
| #endif | | #endif | |
| | | | |
| eT old_avg_log_p = -Datum<eT>::inf; | | eT old_avg_log_p = -Datum<eT>::inf; | |
| | | | |
| for(uword iter=1; iter <= max_iter; ++iter) | | for(uword iter=1; iter <= max_iter; ++iter) | |
| { | | { | |
| init_constants(); | | init_constants(); | |
| | | | |
|
| em_update_params(X, t_boundary, t_acc_means, t_acc_dcovs, t_acc_norm_lh
oods, t_gaus_log_lhoods, t_progress_log_lhood); | | em_update_params(X, boundaries, t_acc_means, t_acc_dcovs, t_acc_norm_lh
oods, t_gaus_log_lhoods, t_progress_log_lhood); | |
| | | | |
| em_fix_params(var_floor); | | em_fix_params(var_floor); | |
| | | | |
| const eT new_avg_log_p = mean(t_progress_log_lhood); | | const eT new_avg_log_p = mean(t_progress_log_lhood); | |
| | | | |
| if(verbose) | | if(verbose) | |
| { | | { | |
| get_stream_err2() << "gmm_diag::learn(): EM: iteration: "; | | get_stream_err2() << "gmm_diag::learn(): EM: iteration: "; | |
| get_stream_err2().unsetf(ios::scientific); | | get_stream_err2().unsetf(ios::scientific); | |
| get_stream_err2().setf(ios::fixed); | | get_stream_err2().setf(ios::fixed); | |
| | | | |
| skipping to change at line 1651 | | skipping to change at line 1790 | |
| | | | |
| return true; | | return true; | |
| } | | } | |
| | | | |
| template<typename eT> | | template<typename eT> | |
| inline | | inline | |
| void | | void | |
| gmm_diag<eT>::em_update_params | | gmm_diag<eT>::em_update_params | |
| ( | | ( | |
| const Mat<eT>& X, | | const Mat<eT>& X, | |
|
| const field< uvec >& t_boundary, | | const umat& boundaries, | |
| field< Mat<eT> >& t_acc_means, | | field< Mat<eT> >& t_acc_means, | |
| field< Mat<eT> >& t_acc_dcovs, | | field< Mat<eT> >& t_acc_dcovs, | |
| field< Col<eT> >& t_acc_norm_lhoods, | | field< Col<eT> >& t_acc_norm_lhoods, | |
| field< Col<eT> >& t_gaus_log_lhoods, | | field< Col<eT> >& t_gaus_log_lhoods, | |
| Col<eT>& t_progress_log_lhood | | Col<eT>& t_progress_log_lhood | |
| ) | | ) | |
| { | | { | |
| arma_extra_debug_sigprint(); | | arma_extra_debug_sigprint(); | |
| | | | |
|
| const uword n_threads = t_boundary.n_elem; | | const uword n_threads = boundaries.n_cols; | |
| | | | |
| // em_generate_acc() is the "map" operation, which produces partial accum
ulators for means, diagonal covariances and hefts | | // em_generate_acc() is the "map" operation, which produces partial accum
ulators for means, diagonal covariances and hefts | |
| | | | |
| #if defined(_OPENMP) | | #if defined(_OPENMP) | |
| { | | { | |
| #pragma omp parallel for | | #pragma omp parallel for | |
| for(uword t=0; t<n_threads; t++) | | for(uword t=0; t<n_threads; t++) | |
| { | | { | |
|
| const uvec& boundary = t_boundary[t]; | | Mat<eT>& acc_means = t_acc_means[t]; | |
| Mat<eT>& acc_means = t_acc_means[t]; | | Mat<eT>& acc_dcovs = t_acc_dcovs[t]; | |
| Mat<eT>& acc_dcovs = t_acc_dcovs[t]; | | Col<eT>& acc_norm_lhoods = t_acc_norm_lhoods[t]; | |
| Col<eT>& acc_norm_lhoods = t_acc_norm_lhoods[t]; | | Col<eT>& gaus_log_lhoods = t_gaus_log_lhoods[t]; | |
| Col<eT>& gaus_log_lhoods = t_gaus_log_lhoods[t]; | | eT& progress_log_lhood = t_progress_log_lhood[t]; | |
| eT& progress_log_lhood = t_progress_log_lhood[t]; | | | |
| | | | |
|
| em_generate_acc(X, boundary, acc_means, acc_dcovs, acc_norm_lhoods, g
aus_log_lhoods, progress_log_lhood); | | em_generate_acc(X, boundaries.at(0,t), boundaries.at(1,t), acc_means,
acc_dcovs, acc_norm_lhoods, gaus_log_lhoods, progress_log_lhood); | |
| } | | } | |
| } | | } | |
| #else | | #else | |
| { | | { | |
|
| em_generate_acc(X, t_boundary[0], t_acc_means[0], t_acc_dcovs[0], t_acc
_norm_lhoods[0], t_gaus_log_lhoods[0], t_progress_log_lhood[0]); | | em_generate_acc(X, boundaries.at(0,0), boundaries.at(1,0), t_acc_means[
0], t_acc_dcovs[0], t_acc_norm_lhoods[0], t_gaus_log_lhoods[0], t_progress_
log_lhood[0]); | |
| } | | } | |
| #endif | | #endif | |
| | | | |
| const uword N_dims = means.n_rows; | | const uword N_dims = means.n_rows; | |
| const uword N_gaus = means.n_cols; | | const uword N_gaus = means.n_cols; | |
| | | | |
| Mat<eT>& final_acc_means = t_acc_means[0]; | | Mat<eT>& final_acc_means = t_acc_means[0]; | |
| Mat<eT>& final_acc_dcovs = t_acc_dcovs[0]; | | Mat<eT>& final_acc_dcovs = t_acc_dcovs[0]; | |
| | | | |
| Col<eT>& final_acc_norm_lhoods = t_acc_norm_lhoods[0]; | | Col<eT>& final_acc_norm_lhoods = t_acc_norm_lhoods[0]; | |
| | | | |
| skipping to change at line 1734 | | skipping to change at line 1872 | |
| } | | } | |
| } | | } | |
| } | | } | |
| | | | |
| template<typename eT> | | template<typename eT> | |
| inline | | inline | |
| void | | void | |
| gmm_diag<eT>::em_generate_acc | | gmm_diag<eT>::em_generate_acc | |
| ( | | ( | |
| const Mat<eT>& X, | | const Mat<eT>& X, | |
|
| const uvec& boundary, | | const uword start_index, | |
| | | const uword end_index, | |
| Mat<eT>& acc_means, | | Mat<eT>& acc_means, | |
| Mat<eT>& acc_dcovs, | | Mat<eT>& acc_dcovs, | |
| Col<eT>& acc_norm_lhoods, | | Col<eT>& acc_norm_lhoods, | |
| Col<eT>& gaus_log_lhoods, | | Col<eT>& gaus_log_lhoods, | |
| eT& progress_log_lhood | | eT& progress_log_lhood | |
| ) | | ) | |
| const | | const | |
| { | | { | |
| arma_extra_debug_sigprint(); | | arma_extra_debug_sigprint(); | |
| | | | |
| | | | |
| skipping to change at line 1759 | | skipping to change at line 1898 | |
| | | | |
| acc_norm_lhoods.zeros(); | | acc_norm_lhoods.zeros(); | |
| gaus_log_lhoods.zeros(); | | gaus_log_lhoods.zeros(); | |
| | | | |
| const uword N_dims = means.n_rows; | | const uword N_dims = means.n_rows; | |
| const uword N_gaus = means.n_cols; | | const uword N_gaus = means.n_cols; | |
| | | | |
| const eT* log_hefts_mem = log_hefts.memptr(); | | const eT* log_hefts_mem = log_hefts.memptr(); | |
| eT* gaus_log_lhoods_mem = gaus_log_lhoods.memptr(); | | eT* gaus_log_lhoods_mem = gaus_log_lhoods.memptr(); | |
| | | | |
|
| const uword start_index = boundary[0]; | | | |
| const uword end_index = boundary[1]; | | | |
| | | | |
| for(uword i=start_index; i <= end_index; i++) | | for(uword i=start_index; i <= end_index; i++) | |
| { | | { | |
| const eT* x = X.colptr(i); | | const eT* x = X.colptr(i); | |
| | | | |
| for(uword g=0; g < N_gaus; ++g) | | for(uword g=0; g < N_gaus; ++g) | |
| { | | { | |
| gaus_log_lhoods_mem[g] = internal_scalar_log_p(x, g) + log_hefts_mem[
g]; | | gaus_log_lhoods_mem[g] = internal_scalar_log_p(x, g) + log_hefts_mem[
g]; | |
| } | | } | |
| | | | |
| eT log_lhood_sum = gaus_log_lhoods_mem[0]; | | eT log_lhood_sum = gaus_log_lhoods_mem[0]; | |
| | | | |
End of changes. 35 change blocks. |
| 102 lines changed or deleted | | 238 lines changed or added | |
|