FORCES
FORtran lib for Comp. Env. Sys.
|
Approximates the probability density function (PDF). More...
Public Member Functions | |
real(dp) function, dimension(:), allocatable | kernel_density_1d_dp (ix, h, silverman, xout, mask, nodata) |
real(sp) function, dimension(:), allocatable | kernel_density_1d_sp (ix, h, silverman, xout, mask, nodata) |
Approximates the probability density function (PDF).
Approximates the probability density function (PDF) to a given 1D data set using a Gaussian kernel.
The bandwith of the kernel can be pre-determined using the function kernel_density_h and specified by the optional argument h. If h is not given the default method to approximate the bandwith h is Silverman''s rule-of-thumb
\[ h = \frac{4}{3}^{0.2} n^{-0.2} \sigma_x \]
where n is the number of given data points and \( \sigma \) is the standard deviation of the data.
If the optional argument silverman is set to false, the cross-validation method described by Scott et al. (2005) is applied. Therefore, the likelihood of a given h is maximized using the Nelder-Mead algorithm nelminrange. For large data sets this might be time consuming and should be performed aforehand using the function kernel_density_h.
The dataset x can be single or double precision. The result will have the same numerical precision.
If the PDF for other datapoints than x is needed the optional argument xout has to be specified. The result will than be of the same size and precision as xout.
Example
-> see also example in test directory
Literature
[in] | real(sp/dp) :: x(:) | \( x_i \) 1D-array with data points |
[in] | real(sp/dp), optional :: h | Bandwith of kernel. If present, argument silverman is ignored. If not present, the bandwith will be approximated first. |
[in] | logical, optional :: silverman | By default Silverman''s Rule-of-thumb will be used to approximate the bandwith of the kernel (silverman=true). If silverman=false the Cross-Validation approach is used to estimate the bandwidth. |
[in] | real(sp/dp), optional :: xout(:) | If present, the PDF will be approximated at this arguments, otherwise the PDF is approximated at x. |
[in] | logical, optional :: mask(:) | mask x values at calculation. if not xout given, then kernel estimates will have nodata value. |
[in] | real(sp/dp), optional :: nodata | if mask and not xout given, then masked data will have nodata kernel estimate. |
Definition at line 206 of file mo_kernel.f90.
real(dp) function, dimension(:), allocatable mo_kernel::kernel_density::kernel_density_1d_dp | ( | real(dp), dimension(:), intent(in) | ix, |
real(dp), intent(in), optional | h, | ||
logical, intent(in), optional | silverman, | ||
real(dp), dimension(:), intent(in), optional | xout, | ||
logical, dimension(:), intent(in), optional | mask, | ||
real(dp), intent(in), optional | nodata | ||
) |
Definition at line 861 of file mo_kernel.f90.
real(sp) function, dimension(:), allocatable mo_kernel::kernel_density::kernel_density_1d_sp | ( | real(sp), dimension(:), intent(in) | ix, |
real(sp), intent(in), optional | h, | ||
logical, intent(in), optional | silverman, | ||
real(sp), dimension(:), intent(in), optional | xout, | ||
logical, dimension(:), intent(in), optional | mask, | ||
real(sp), intent(in), optional | nodata | ||
) |
Definition at line 957 of file mo_kernel.f90.