This is an implementation of the probability density function and a random number generator for Student's t-distribution. (See the Wikipedia article [http://en.wikipedia.org/wiki/Student%27s_t-distribution].)

The function cdf-students-t is already implemented in [tcllib], as is a procedure for calculating the t-statistic by calculating the inverse CDF.

This code implements pdf-students-t and random-students-t. Included is the code for the functions as well as some test code. The test code uses [Simple Test].

**The Library (tdist.tcl)**
======
    package require math
    namespace import math::ln_Gamma

    source simpletest.tcl

    namespace eval tdist {
        variable pi [expr {atan2(0, -1)}]

        st::addproc pdf-students-t {
            degrees     ;# Degrees of freedom (allow non-integer)
            value       ;# svalue where dist is to be evaluated
        } where {
            {{1 0.1} ~ 0.315158303152268}
            {{0.5 0.1} ~ 0.265700672177405}
            {{4 3.2} ~ 0.0156821741652879}
            {{3 2.0} ~ 0.0675096606638929}
            {{3 7.5} ~ 0.000942291548015668}
        }
        
        st::addproc random-students-t {
            degrees ;# Degrees of freedom: may be non-integer, but >= 1
            number  ;# Number of random deviates to produce
        } where {} ;# No tests: compare visually as histogram
    }

    proc tdist::pdf-students-t {degrees value} {
        variable pi
        
        set nplus1over2 [expr {0.5 * ($degrees + 1)}]
        set f1 [expr {exp([ln_Gamma $nplus1over2] - \
            [ln_Gamma [expr {$nplus1over2 - 0.5}]])}]
        set f2 [expr {1.0/sqrt($degrees * $pi)}]
        
        expr {$f1 * $f2 * pow(1.0 + $value * $value/double($degrees), -$nplus1over2)} 
    }

    # Use method from Appendix 4.3 in Dagpunar, J.S.,
    # "Simulation and Monte Carlo: With Applications in
    # Finance and MCMC"
    proc tdist::random-students-t {degrees number} {
        variable pi
        
        if {$degrees < 1} {
            error "random-students-t: Degrees of freedom must be >= 1"
        }
        
        set dd [expr {double($degrees)}]
        set k [expr {2.0/($dd - 1.0)}]
        
        for {set i 0} {$i < $number} {incr i} {
            set r1 [expr {rand()}]
            if {$degrees > 1} {
                set r2 [expr {rand()}]
                set c [expr {cos(2.0 * $pi * $r2)}]
                lappend retval [expr {sqrt($dd/ \
                    (1.0/(1.0 - pow($r1, $k)) \
                    - $c * $c)) * $c}]
            } else {
                lappend retval [expr {tan(0.5 * $pi * ($r1 + $r1 - 1))}]
            }
        }
        set retval
    }
======

**Test Suite (tdist_test.tcl)**
======
    package require math::statistics
    source tdist.tcl

    st::tolerance 1.0e-9
    console show
    namespace eval tdist {
        st::print
    }

    ##########################################################################################
    ##
    ## TESTING
    ##
    ## For random numbers, generate histograms:
    ##
    ##########################################################################################
    canvas .c
    pack .c -side top
    frame .f
    pack .f -side bottom
    label .f.dofl -text "dof"
    entry .f.dofe -textvariable dof
    pack .f.dofl -side left
    pack .f.dofe -side left
    button .f.run -text "Run" -command runtest
    pack .f.run -side left

    proc runtest {} {
        set step 0.5
        
        set vals [tdist::random-students-t $::dof 5000]
        set remainder 5000
        for {set x -10.0} {$x < 10.0} {set x [expr {$x + $step}]} {
            lappend bins $x
            set distval [tdist::pdf-students-t $::dof [expr {$x - 0.5 * $step}]]
            # Total trials (5000) * step * value
            set distval [expr {int(5000 * $step * $distval)}]
            lappend distcounts $distval
        }
        # Assume none are left
        lappend distcounts 0.0
        set bincounts [::math::statistics::histogram $bins $vals]
        .c delete hist
        .c delete dist
        math::statistics::plot-scale .c -7 7 0 [math::statistics::max $bincounts]
        math::statistics::plot-histogram .c $bincounts $bins hist
        math::statistics::plot-histogram .c $distcounts $bins dist
        .c itemconfigure dist -fill {} -outline green
    }

    set dof 3
    runtest
======


<<categories>> Mathematics | Statistics