From 67c005341d207cc66d3bcef6f3fbb2e3108b4234 Mon Sep 17 00:00:00 2001 From: "per.andreas.brodtkorb" Date: Wed, 26 Jan 2011 15:16:14 +0000 Subject: [PATCH] Fixed a bug in nlogps when ties occuring in the data. --- pywafo/src/wafo/source/mreg/quadrmod.mod | 2 +- pywafo/src/wafo/source/mreg/svd.mod | 40 ++++++++++++------------ pywafo/src/wafo/stats/distributions.py | 4 +-- pywafo/src/wafo/stats/estimation.py | 21 +++++++++++-- 4 files changed, 42 insertions(+), 25 deletions(-) diff --git a/pywafo/src/wafo/source/mreg/quadrmod.mod b/pywafo/src/wafo/source/mreg/quadrmod.mod index 5f4614f..2cb84fd 100644 --- a/pywafo/src/wafo/source/mreg/quadrmod.mod +++ b/pywafo/src/wafo/source/mreg/quadrmod.mod @@ -1,4 +1,4 @@ -GFORTRAN module version '0' created from mregmodule.f on Wed Aug 05 19:15:05 2009 +GFORTRAN module version '0' created from mregmodule.f on Tue Jan 25 23:52:18 2011 MD5:447301769c212f228b6cfa086ba1d48a -- If you edit this, you'll get what you deserve. (() () () () () () () () () () () () () () () () () () () () () () () () diff --git a/pywafo/src/wafo/source/mreg/svd.mod b/pywafo/src/wafo/source/mreg/svd.mod index 7abb5f6..5fb3a0e 100644 --- a/pywafo/src/wafo/source/mreg/svd.mod +++ b/pywafo/src/wafo/source/mreg/svd.mod @@ -1,5 +1,5 @@ -GFORTRAN module version '0' created from dsvdc.f on Wed Aug 05 19:15:05 2009 -MD5:06f4ab5dbb3a45df847b6d37183b2196 -- If you edit this, you'll get what you deserve. +GFORTRAN module version '0' created from dsvdc.f on Tue Jan 25 20:12:18 2011 +MD5:504db75f3667c360354623f37288fa05 -- If you edit this, you'll get what you deserve. (() () () () () () () () () () () () () () () () () () () () () () () () () () ()) @@ -32,19 +32,22 @@ PROCEDURE UNKNOWN-INTENT UNKNOWN-PROC UNKNOWN UNKNOWN FUNCTION) ( UNKNOWN 0 0 0 UNKNOWN ()) 0 0 () () 32 () () () 0 0) 33 'svd' 'svd' 'svd' 1 ((MODULE UNKNOWN-INTENT UNKNOWN-PROC UNKNOWN UNKNOWN) (UNKNOWN 0 0 0 UNKNOWN ()) 0 0 () () 0 () () () 0 0) -15 'ds' '' 'ds' 11 ((VARIABLE OUT UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) ( -REAL 8 0 0 REAL ()) 0 0 () () 0 () () () 0 0) 12 'da' '' 'da' 11 ((VARIABLE INOUT UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) (REAL 8 0 0 REAL ()) 0 0 () () 0 () () () 0 0) 13 'db' '' 'db' 11 ((VARIABLE INOUT UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) (REAL 8 0 0 REAL ()) 0 0 () () 0 () () () 0 0) 14 'dc' '' 'dc' 11 ((VARIABLE OUT UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) ( REAL 8 0 0 REAL ()) 0 0 () () 0 () () () 0 0) +15 'ds' '' 'ds' 11 ((VARIABLE OUT UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) ( +REAL 8 0 0 REAL ()) 0 0 () () 0 () () () 0 0) 29 'n' '' 'n' 28 ((VARIABLE IN UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) ( INTEGER 4 0 0 INTEGER ()) 0 0 () () 0 () () () 0 0) 30 'dx' '' 'dx' 28 ((VARIABLE INOUT UNKNOWN-PROC UNKNOWN UNKNOWN DIMENSION DUMMY) (REAL 8 0 0 REAL ()) 0 0 () (1 ASSUMED_SIZE (CONSTANT ( INTEGER 4 0 0 INTEGER ()) 0 '1') ()) 0 () () () 0 0) +31 'dy' '' 'dy' 28 ((VARIABLE INOUT UNKNOWN-PROC UNKNOWN UNKNOWN +DIMENSION DUMMY) (REAL 8 0 0 REAL ()) 0 0 () (1 ASSUMED_SIZE (CONSTANT ( +INTEGER 4 0 0 INTEGER ()) 0 '1') ()) 0 () () () 0 0) 5 'n' '' 'n' 4 ((VARIABLE IN UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) ( INTEGER 4 0 0 INTEGER ()) 0 0 () () 0 () () () 0 0) 6 'dx' '' 'dx' 4 ((VARIABLE INOUT UNKNOWN-PROC UNKNOWN UNKNOWN DIMENSION @@ -53,6 +56,18 @@ DUMMY) (REAL 8 0 0 REAL ()) 0 0 () (1 ASSUMED_SIZE (CONSTANT (INTEGER 4 7 'dy' '' 'dy' 4 ((VARIABLE INOUT UNKNOWN-PROC UNKNOWN UNKNOWN DIMENSION DUMMY) (REAL 8 0 0 REAL ()) 0 0 () (1 ASSUMED_SIZE (CONSTANT (INTEGER 4 0 0 INTEGER ()) 0 '1') ()) 0 () () () 0 0) +8 'c' '' 'c' 4 ((VARIABLE IN UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) (REAL 8 +0 0 REAL ()) 0 0 () () 0 () () () 0 0) +9 's' '' 's' 4 ((VARIABLE IN UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) (REAL 8 +0 0 REAL ()) 0 0 () () 0 () () () 0 0) +18 'x' '' 'x' 17 ((VARIABLE INOUT UNKNOWN-PROC UNKNOWN UNKNOWN DIMENSION +DUMMY) (REAL 8 0 0 REAL ()) 0 0 () (2 ASSUMED_SHAPE (CONSTANT (INTEGER 4 +0 0 INTEGER ()) 0 '1') () (CONSTANT (INTEGER 4 0 0 INTEGER ()) 0 '1') ()) +0 () () () 0 0) +19 'n' '' 'n' 17 ((VARIABLE IN UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) ( +INTEGER 4 0 0 INTEGER ()) 0 0 () () 0 () () () 0 0) +20 'p' '' 'p' 17 ((VARIABLE IN UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) ( +INTEGER 4 0 0 INTEGER ()) 0 0 () () 0 () () () 0 0) 21 's' '' 's' 17 ((VARIABLE OUT UNKNOWN-PROC UNKNOWN UNKNOWN DIMENSION DUMMY) (REAL 8 0 0 REAL ()) 0 0 () (1 ASSUMED_SHAPE (CONSTANT (INTEGER 4 0 0 INTEGER ()) 0 '1') ()) 0 () () () 0 0) @@ -63,29 +78,14 @@ DUMMY) (REAL 8 0 0 REAL ()) 0 0 () (1 ASSUMED_SHAPE (CONSTANT (INTEGER 4 DUMMY) (REAL 8 0 0 REAL ()) 0 0 () (2 ASSUMED_SHAPE (CONSTANT (INTEGER 4 0 0 INTEGER ()) 0 '1') () (CONSTANT (INTEGER 4 0 0 INTEGER ()) 0 '1') ()) 0 () () () 0 0) -25 'job' '' 'job' 17 ((VARIABLE IN UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) ( -INTEGER 4 0 0 INTEGER ()) 0 0 () () 0 () () () 0 0) 24 'v' '' 'v' 17 ((VARIABLE OUT UNKNOWN-PROC UNKNOWN UNKNOWN DIMENSION DUMMY) (REAL 8 0 0 REAL ()) 0 0 () (2 ASSUMED_SHAPE (CONSTANT (INTEGER 4 0 0 INTEGER ()) 0 '1') () (CONSTANT (INTEGER 4 0 0 INTEGER ()) 0 '1') ()) 0 () () () 0 0) -18 'x' '' 'x' 17 ((VARIABLE INOUT UNKNOWN-PROC UNKNOWN UNKNOWN DIMENSION -DUMMY) (REAL 8 0 0 REAL ()) 0 0 () (2 ASSUMED_SHAPE (CONSTANT (INTEGER 4 -0 0 INTEGER ()) 0 '1') () (CONSTANT (INTEGER 4 0 0 INTEGER ()) 0 '1') ()) -0 () () () 0 0) -19 'n' '' 'n' 17 ((VARIABLE IN UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) ( -INTEGER 4 0 0 INTEGER ()) 0 0 () () 0 () () () 0 0) -20 'p' '' 'p' 17 ((VARIABLE IN UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) ( +25 'job' '' 'job' 17 ((VARIABLE IN UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) ( INTEGER 4 0 0 INTEGER ()) 0 0 () () 0 () () () 0 0) 26 'info' '' 'info' 17 ((VARIABLE OUT UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) (INTEGER 4 0 0 INTEGER ()) 0 0 () () 0 () () () 0 0) -31 'dy' '' 'dy' 28 ((VARIABLE INOUT UNKNOWN-PROC UNKNOWN UNKNOWN -DIMENSION DUMMY) (REAL 8 0 0 REAL ()) 0 0 () (1 ASSUMED_SIZE (CONSTANT ( -INTEGER 4 0 0 INTEGER ()) 0 '1') ()) 0 () () () 0 0) -8 'c' '' 'c' 4 ((VARIABLE IN UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) (REAL 8 -0 0 REAL ()) 0 0 () () 0 () () () 0 0) -9 's' '' 's' 4 ((VARIABLE IN UNKNOWN-PROC UNKNOWN UNKNOWN DUMMY) (REAL 8 -0 0 REAL ()) 0 0 () () 0 () () () 0 0) ) ('dp' 0 2 'drot1' 0 3 'drotg' 0 10 'dsvdc' 0 16 'dswap1' 0 27 diff --git a/pywafo/src/wafo/stats/distributions.py b/pywafo/src/wafo/stats/distributions.py index dcd3f7a..5f8f74b 100644 --- a/pywafo/src/wafo/stats/distributions.py +++ b/pywafo/src/wafo/stats/distributions.py @@ -2032,7 +2032,7 @@ class rv_continuous(rv_generic): if (any(cond0)): return inf else: - linfo = numpy.finfo(float) + #linfo = numpy.finfo(float) realmax = floatinfo.machar.xmax lowertail = True @@ -2059,7 +2059,7 @@ class rv_continuous(rv_generic): i_tie = nonzero(tie) tiedata = x[i_tie] - logD[i_tie + 1] = log(self._pdf(tiedata, *args)) + log(scale) + logD[i_tie + 1] = log(self._pdf(tiedata, *args)) - log(scale) finiteD = numpy.isfinite(logD) nonfiniteD = 1 - finiteD diff --git a/pywafo/src/wafo/stats/estimation.py b/pywafo/src/wafo/stats/estimation.py index 7172183..4b4988d 100644 --- a/pywafo/src/wafo/stats/estimation.py +++ b/pywafo/src/wafo/stats/estimation.py @@ -901,8 +901,25 @@ class FitDistribution(rv_frozen): If so the histogram should resemble the model density. Other distribution types will introduce deviations in the plot. ''' - - bin, limits = numpy.histogram(self.data, normed=True) #, new=True) + odd = False + x = np.atleast_1d(self.data) + n = np.ceil(4 * np.sqrt(np.sqrt(len(x)))) + + mn = x.min() + mx = x.max() + d = (mx - mn) / n * 2 + e = np.floor(np.log(d) / np.log(10)); + m = np.floor(d / 10 ** e) + if m > 5: + m = 5 + elif m > 2: + m = 2 + + d = m * 10 ** e + mn = (np.floor(mn / d) - 1) * d - odd * d / 2 + mx = (np.ceil(mx / d) + 1) * d + odd * d / 2 + limits = np.arange(mn, mx, d) + bin, limits = numpy.histogram(self.data, bins=limits, normed=True) #, new=True) limits.shape = (-1, 1) xx = limits.repeat(3, axis=1) xx.shape = (-1,)