function dirichlet_mix_sample_test ( ) %*****************************************************************************80 % %% DIRICHLET_MIX_SAMPLE_TEST tests DIRICHLET_MIX_SAMPLE. % % Licensing: % % This code is distributed under the GNU LGPL license. % % Modified: % % 13 April 2016 % % Author: % % John Burkardt % comp_num = 2; elem_num = 3; sample_num = 1000; seed = 123456789; fprintf ( 1, '\n' ); fprintf ( 1, 'DIRICHLET_MIX_SAMPLE_TEST\n' ); fprintf ( 1, ' DIRICHLET_MIX_SAMPLE samples the Dirichlet Mix distribution;\n' ); fprintf ( 1, ' DIRICHLET_MIX_MEAN computes the Dirichlet Mix mean;\n' ); a(1:3,1:2) = [ 0.250, 1.500; ... 0.500, 0.500; ... 1.250, 2.000 ]; comp_weight(1:2) = [ 1.0, 2.0 ]; check = dirichlet_mix_check ( comp_num, elem_num, a, comp_weight ); if ( ~ check ); fprintf ( 1, '\n' ); fprintf ( 1, 'DIRICHLET_MIX_SAMPLE_TEST - Fatal error!\n' ); fprintf ( 1, ' The parameters are not legal.\n' ); return end fprintf ( 1, '\n' ); fprintf ( 1, ' Number of elements ELEM_NUM = %6d\n', elem_num ); fprintf ( 1, ' Number of components COMP_NUM = %6d\n', comp_num ); r8mat_print ( elem_num, comp_num, a, ' PDF parameters A(ELEM,COMP):' ); r8vec_print ( comp_num, comp_weight, ' Component weights:' ); mean(1:elem_num) = dirichlet_mix_mean ( comp_num, elem_num, a, comp_weight ); r8vec_print ( elem_num, mean, ' PDF mean:' ); for sample_i = 1 : sample_num [ x(1:elem_num,sample_i), comp, seed ] = dirichlet_mix_sample ( ... comp_num, elem_num, a, comp_weight, seed ); end xmax(1:elem_num) = r8row_max ( elem_num, sample_num, x ); xmin(1:elem_num) = r8row_min ( elem_num, sample_num, x ); mean(1:elem_num) = r8row_mean ( elem_num, sample_num, x ); variance(1:elem_num) = r8row_variance ( elem_num, sample_num, x ); fprintf ( 1, '\n' ); fprintf ( 1, ' Sample size = %6d\n', sample_num ); fprintf ( 1, '\n' ); fprintf ( 1, ' Observed Min, Max, Mean, Variance:\n' ); fprintf ( 1, '\n' ); for elem_i = 1 : elem_num fprintf ( 1, ' %6d %14g %14g %14g %14g\n', elem_i, xmin(elem_i), ... xmax(elem_i), mean(elem_i), variance(elem_i) ); end return end