#!/usr/local/bin/webwork-perl =head1 Statistics Macros =head3 Normal distribution =pod Usage: normal(a, b, mean=>0, deviation=>1); Computes the probability of x being in the interval (a,b) for normal distribution. The first two arguments are required. Use '-infty' for negative infinity, and 'infty' or '+infty' for positive infinity. The mean and deviation are optional, and are 0 and 1 respectively by default. =cut sub normal { my \$a = shift; my \$b = shift; my %options = @_; my \$mean = \$options{'mean'} if defined (\$options{'mean'}); \$mean = 0 unless defined \$mean; my \$deviation = \$options{'deviation'} if defined (\$options{'deviation'}); \$deviation = 1 unless defined \$deviation; if (\$deviation <= 0) { warn 'Deviation must be a positive number.'; return; } if ( \$a eq '-infty' ) { \$a = -6*\$deviation + \$mean; } if ( \$b eq 'infty' ) { \$b = 6*\$deviation + \$mean; } if ( \$b eq '+infty' ) { \$b = 6*\$deviation + \$mean; } my \$z_score_of_a = (\$a - \$mean)/\$deviation; my \$z_score_of_b = (\$b - \$mean)/\$deviation; my \$function = sub { my \$x=shift; \$E**(-\$x**2/2)/sqrt(2*\$PI); }; my \$prob = romberg(\$function, \$z_score_of_a, \$z_score_of_b); \$prob; } =head3 "Inverse" of normal distribution =pod Usage: inv_normal(prob, mean=>0, deviation=>1); Computes the positive number b such that the probability of x being in the interval (0,b) is equal to the given probability (first argument). The mean and deviation are optional, and are 0 and 1 respectively by default. Caution: since students may use tables, they may only be able to provide the answer correct to 2 or 3 decimal places. Use tolerance when evaluating answers. =cut sub inv_normal { my \$prob = shift; my %options = @_; my \$mean = \$options{'mean'} if defined (\$options{'mean'}); \$mean = 0 unless defined \$mean; my \$deviation = \$options{'deviation'} if defined (\$options{'deviation'}); \$deviation = 1 unless defined \$deviation; if (\$deviation <= 0) { warn 'Deviation must be a positive number.'; return; } my \$function = sub { my \$x=shift; \$E**(-\$x**2/2)/sqrt(2*\$PI); }; my \$z_score_of_b = inv_romberg(\$function, 0, \$prob); my \$b = \$z_score_of_b * \$deviation + \$mean; \$b; } ########################################## 1;