Python scipy.stats.norm() Examples

def test_qqplot(self):
        """Test qqplot()"""
        np.random.seed(123)
        x = np.random.normal(size=50)
        x_ln = np.random.lognormal(size=50)
        x_exp = np.random.exponential(size=50)
        ax = qqplot(x, dist='norm')
        assert isinstance(ax, matplotlib.axes.Axes)
        _, (ax1, ax2) = plt.subplots(1, 2, figsize=(9, 4))
        qqplot(x_exp, dist='expon', ax=ax2)
        mean, std = 0, 0.8
        qqplot(x, dist=stats.norm, sparams=(mean, std), confidence=False)
        # For lognormal distribution, the shape parameter must be specified
        ax = qqplot(x_ln, dist='lognorm', sparams=(1))
        assert isinstance(ax, matplotlib.axes.Axes)
        # Error: required parameters are not specified
        with pytest.raises(ValueError):
            qqplot(x_ln, dist='lognorm', sparams=())
        plt.close('all') 

def setUp_configure(self):
        from scipy import stats
        self.dist = distributions.Normal
        self.scipy_dist = stats.norm

        self.test_targets = set([
            'batch_shape', 'cdf', 'entropy', 'event_shape', 'icdf', 'log_cdf',
            'log_prob', 'log_survival', 'mean', 'prob', 'sample', 'stddev',
            'support', 'survival', 'variance'])

        loc = utils.force_array(
            numpy.random.uniform(-1, 1, self.shape).astype(numpy.float32))
        if self.log_scale_option:
            log_scale = utils.force_array(
                numpy.random.uniform(-1, 1, self.shape).astype(numpy.float32))
            scale = numpy.exp(log_scale)
            self.params = {'loc': loc, 'log_scale': log_scale}
            self.scipy_params = {'loc': loc, 'scale': scale}
        else:
            scale = utils.force_array(numpy.exp(
                numpy.random.uniform(-1, 1, self.shape)).astype(numpy.float32))
            self.params = {'loc': loc, 'scale': scale}
            self.scipy_params = {'loc': loc, 'scale': scale} 

def __init__(self,
                 nbases,
                 Xdim,
                 mean=Parameter(norm_dist(), Bound()),
                 lenscale=Parameter(gamma(1.), Positive()),
                 regularizer=None,
                 random_state=None
                 ):
        """See this class's docstring."""
        self.random_state = random_state  # for repr
        self._random = check_random_state(random_state)
        self._init_dims(nbases, Xdim)
        self._params = [self._init_param(mean),
                        self._init_param(lenscale)]
        self._init_matrices()
        super(_LengthScaleBasis, self).__init__(regularizer) 

def test_mixture_rvs_fixed(self):
        mix = MixtureDistribution()
        np.random.seed(1234)
        res = mix.rvs([.15,.85], 50, dist=[stats.norm, stats.norm], kwargs =
                (dict(loc=1,scale=.5),dict(loc=-1,scale=.5)))
        npt.assert_almost_equal(
                res,
                np.array([-0.5794956 , -1.72290504, -1.70098664, -1.0504591 ,
                            -1.27412122,-1.07230975, -0.82298983, -1.01775651,
                            -0.71713085,-0.2271706 ,-1.48711817, -1.03517244,
                            -0.84601557, -1.10424938, -0.48309963,-2.20022682,
                            0.01530181,  1.1238961 , -1.57131564, -0.89405831,
                            -0.64763969, -1.39271761,  0.55142161, -0.76897013,
                            -0.64788589,-0.73824602, -1.46312716,  0.00392148,
                            -0.88651873, -1.57632955,-0.68401028, -0.98024366,
                            -0.76780384,  0.93160258,-2.78175833,-0.33944719,
                            -0.92368472, -0.91773523, -1.21504785, -0.61631563,
                            1.0091446 , -0.50754008,  1.37770699, -0.86458208,
                            -0.3040069 ,-0.96007884,  1.10763429, -1.19998229,
                            -1.51392528, -1.29235911])) 

def test_cdf_ppf(self):
        values = np.array([0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5,
                           5.0, 5.5, 6.0, 6.5, 7.0, 7.5, 8.0, 8.5, 9.0, 9.5])
        cdf_values = np.asarray([0.0/25.0, 0.0/25.0, 0.0/25.0, 0.5/25.0,
                                 1.0/25.0, 2.0/25.0, 3.0/25.0, 4.5/25.0,
                                 6.0/25.0, 8.0/25.0, 10.0/25.0, 12.5/25.0,
                                 15.0/25.0, 17.0/25.0, 19.0/25.0, 20.5/25.0,
                                 22.0/25.0, 23.5/25.0, 25.0/25.0, 25.0/25.0])
        assert_allclose(self.template.cdf(values), cdf_values)
        # First three and last two values in cdf_value are not unique
        assert_allclose(self.template.ppf(cdf_values[2:-1]), values[2:-1])

        # Test of cdf and ppf are inverse functions
        x = np.linspace(1.0, 9.0, 100)
        assert_allclose(self.template.ppf(self.template.cdf(x)), x)
        x = np.linspace(0.0, 1.0, 100)
        assert_allclose(self.template.cdf(self.template.ppf(x)), x)

        x = np.linspace(-2, 2, 10)
        assert_allclose(self.norm_template.cdf(x),
                        stats.norm.cdf(x, loc=1.0, scale=2.5), rtol=0.1) 

def test_pdf(self):
        values = np.array([0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5,
                           5.0, 5.5, 6.0, 6.5, 7.0, 7.5, 8.0, 8.5, 9.0, 9.5])
        pdf_values = np.asarray([0.0/25.0, 0.0/25.0, 1.0/25.0, 1.0/25.0,
                                 2.0/25.0, 2.0/25.0, 3.0/25.0, 3.0/25.0,
                                 4.0/25.0, 4.0/25.0, 5.0/25.0, 5.0/25.0,
                                 4.0/25.0, 4.0/25.0, 3.0/25.0, 3.0/25.0,
                                 3.0/25.0, 3.0/25.0, 0.0/25.0, 0.0/25.0])
        assert_allclose(self.template.pdf(values), pdf_values)

        # Test explicitly the corner cases:
        # As stated above the pdf in the bin [8,9) is greater than
        # one would naively expect because np.histogram putted the 9
        # into the [8,9) bin.
        assert_almost_equal(self.template.pdf(8.0), 3.0/25.0)
        assert_almost_equal(self.template.pdf(8.5), 3.0/25.0)
        # 9 is outside our defined bins [8,9) hence the pdf is already 0
        # for a continuous distribution this is fine, because a single value
        # does not have a finite probability!
        assert_almost_equal(self.template.pdf(9.0), 0.0/25.0)
        assert_almost_equal(self.template.pdf(10.0), 0.0/25.0)

        x = np.linspace(-2, 2, 10)
        assert_allclose(self.norm_template.pdf(x),
                        stats.norm.pdf(x, loc=1.0, scale=2.5), rtol=0.1) 

def test_compare(self):
        xx = self.res1.support
        kde_vals = [self.res1.evaluate(xi) for xi in xx]
        kde_vals = np.squeeze(kde_vals)  #kde_vals is a "column_list"
        mask_valid = np.isfinite(kde_vals)
        # TODO: nans at the boundaries
        kde_vals[~mask_valid] = 0
        npt.assert_almost_equal(self.res1.density, kde_vals,
                                self.decimal_density)

        # regression test, not compared to another package
        nobs = len(self.res1.endog)
        kern = self.res1.kernel
        v = kern.density_var(kde_vals, nobs)
        v_direct = kde_vals * kern.L2Norm / kern.h / nobs
        npt.assert_allclose(v, v_direct, rtol=1e-10)

        ci = kern.density_confint(kde_vals, nobs)
        crit = 1.9599639845400545 #stats.norm.isf(0.05 / 2)
        hw = kde_vals - ci[:, 0]
        npt.assert_allclose(hw, crit * np.sqrt(v), rtol=1e-10)
        hw = ci[:, 1] - kde_vals
        npt.assert_allclose(hw, crit * np.sqrt(v), rtol=1e-10) 

def test_compute_perfect_model_da1d_not_nan_crpss_quadratic_kwargs(
    PM_da_initialized_1d, PM_da_control_1d
):
    """
    Checks that there are no NaNs on perfect model metrics of 1D time series.
    """
    actual = (
        compute_perfect_model(
            PM_da_initialized_1d.isel(lead=[0]),
            PM_da_control_1d,
            comparison='m2c',
            metric='crpss',
            gaussian=False,
            dim='member',
            tol=1e-6,
            xmin=None,
            xmax=None,
            cdf_or_dist=norm,
        )
        .isnull()
        .any()
    )
    assert not actual 

def __init__(self, predicted_mean, var_pred_mean, var_resid,
                 df=None, dist=None, row_labels=None):
        self.predicted_mean = predicted_mean
        self.var_pred_mean = var_pred_mean
        self.df = df
        self.var_resid = var_resid
        self.row_labels = row_labels

        if dist is None or dist == 'norm':
            self.dist = stats.norm
            self.dist_args = ()
        elif dist == 't':
            self.dist = stats.t
            self.dist_args = (self.df,)
        else:
            self.dist = dist
            self.dist_args = () 

def test_norm_logcdf():
    # Test precision of the logcdf of the normal distribution.
    # This precision was enhanced in ticket 1614.
    x = -np.asarray(list(range(0, 120, 4)))
    # Values from R
    expected = [-0.69314718, -10.36010149, -35.01343716, -75.41067300,
                -131.69539607, -203.91715537, -292.09872100, -396.25241451,
                -516.38564863, -652.50322759, -804.60844201, -972.70364403,
                -1156.79057310, -1356.87055173, -1572.94460885, -1805.01356068,
                -2053.07806561, -2317.13866238, -2597.19579746, -2893.24984493,
                -3205.30112136, -3533.34989701, -3877.39640444, -4237.44084522,
                -4613.48339520, -5005.52420869, -5413.56342187, -5837.60115548,
                -6277.63751711, -6733.67260303]

    assert_allclose(stats.norm().logcdf(x), expected, atol=1e-8)

    # also test the complex-valued code path
    assert_allclose(stats.norm().logcdf(x + 1e-14j).real, expected, atol=1e-8)

    # test the accuracy: d(logcdf)/dx = pdf / cdf \equiv exp(logpdf - logcdf)
    deriv = (stats.norm.logcdf(x + 1e-10j)/1e-10).imag
    deriv_expected = np.exp(stats.norm.logpdf(x) - stats.norm.logcdf(x))
    assert_allclose(deriv, deriv_expected, atol=1e-10) 

def conf_int(self, alpha=.05):
        """
        Returns the confidence intervals of the marginal effects

        Parameters
        ----------
        alpha : float
            Number between 0 and 1. The confidence intervals have the
            probability 1-alpha.

        Returns
        -------
        conf_int : ndarray
            An array with lower, upper confidence intervals for the marginal
            effects.
        """
        _check_at_is_all(self.margeff_options)
        me_se = self.margeff_se
        q = stats.norm.ppf(1 - alpha / 2)
        lower = self.margeff - q * me_se
        upper = self.margeff + q * me_se
        return np.asarray(lzip(lower, upper)) 

def test_star_args(self):
        # test _pdf with only starargs
        # NB: **kwargs of pdf will never reach _pdf
        class _dist_gen(stats.rv_continuous):
            def _pdf(self, x, *args):
                extra_kwarg = args[0]
                return stats.norm._pdf(x) * extra_kwarg

        dist = _dist_gen(shapes='extra_kwarg')
        assert_equal(dist.pdf(0.5, extra_kwarg=33), stats.norm.pdf(0.5)*33)
        assert_equal(dist.pdf(0.5, 33), stats.norm.pdf(0.5)*33)
        assert_raises(TypeError, dist.pdf, 0.5, **dict(xxx=33)) 

def bell_shaped_effect(self, amplitude, width, lag=0, cut=0):
        """Returns coefficients corresponding to a bell shaped relative
        incidence of max value equal to `amplitude`. If `cut` is greater than 0,
        the relative incidence will be null on [`cut`, `n_intervals`]. The
        effect starts at `lag` interval, and lasts `width` intervals.
        """
        self._check_params(lag, width, amplitude, cut)
        if width % 2 == 0:
            width += 1
        effect = norm(0, width / 5).pdf(np.arange(width) - int(width / 2))
        return self._create_risk_curve(effect, amplitude, cut, width, lag) 

def setup_method(self):
        np.random.seed(1234)

        # We have 8 bins
        # [1,2), [2,3), [3,4), [4,5), [5,6), [6,7), [7,8), [8,9)
        # But actually np.histogram will put the last 9 also in the [8,9) bin!
        # Therefore there is a slight difference below for the last bin, from
        # what you might have expected.
        histogram = np.histogram([1, 2, 2, 3, 3, 3, 4, 4, 4, 4, 5, 5, 5, 5, 5,
                                  6, 6, 6, 6, 7, 7, 7, 8, 8, 9], bins=8)
        self.template = stats.rv_histogram(histogram)

        data = stats.norm.rvs(loc=1.0, scale=2.5,size=10000, random_state=123)
        norm_histogram = np.histogram(data, bins=50)
        self.norm_template = stats.rv_histogram(norm_histogram) 

def test_probplot_qq_dist(plot_data):
    fig, ax = plt.subplots()
    norm = stats.norm(*stats.norm.fit(plot_data))
    fig = viz.probplot(plot_data, ax=ax, plottype='qq', dist=norm,
                       datalabel='Test label')
    return fig 

def test_star_args_2(self):
        # test _pdf with named & starargs
        # NB: **kwargs of pdf will never reach _pdf
        class _dist_gen(stats.rv_continuous):
            def _pdf(self, x, offset, *args):
                extra_kwarg = args[0]
                return stats.norm._pdf(x) * extra_kwarg + offset

        dist = _dist_gen(shapes='offset, extra_kwarg')
        assert_equal(dist.pdf(0.5, offset=111, extra_kwarg=33),
                     stats.norm.pdf(0.5)*33 + 111)
        assert_equal(dist.pdf(0.5, 111, 33),
                     stats.norm.pdf(0.5)*33 + 111) 

def test_extra_kwarg(self):
        # **kwargs to _pdf are ignored.
        # this is a limitation of the framework (_pdf(x, *goodargs))
        class _distr_gen(stats.rv_continuous):
            def _pdf(self, x, *args, **kwargs):
                # _pdf should handle *args, **kwargs itself.  Here "handling"
                # is ignoring *args and looking for ``extra_kwarg`` and using
                # that.
                extra_kwarg = kwargs.pop('extra_kwarg', 1)
                return stats.norm._pdf(x) * extra_kwarg

        dist = _distr_gen(shapes='extra_kwarg')
        assert_equal(dist.pdf(1, extra_kwarg=3), stats.norm.pdf(1)) 

def shapes_empty_string(self):
        # shapes='' is equivalent to shapes=None
        class _dist_gen(stats.rv_continuous):
            def _pdf(self, x):
                return stats.norm.pdf(x)

        dist = _dist_gen(shapes='')
        assert_equal(dist.pdf(0.5), stats.norm.pdf(0.5)) 

def test_entropy(self):
        assert_allclose(self.norm_template.entropy(),
                        stats.norm.entropy(loc=1.0, scale=2.5), rtol=0.05) 

def test_logcdf(self):
        # Regression test for gh-5940: sf et al would underflow too early
        x2, mu, sigma = 201.68, 195, 0.149
        assert_allclose(stats.lognorm.sf(x2-mu, s=sigma),
                        stats.norm.sf(np.log(x2-mu)/sigma))
        assert_allclose(stats.lognorm.logsf(x2-mu, s=sigma),
                        stats.norm.logsf(np.log(x2-mu)/sigma)) 

def setUp_configure(self):
        from scipy import stats
        self.dist = lambda **params: distributions.Independent(
            distributions.Normal(**params), self.reinterpreted_batch_ndims)

        self.test_targets = set([
            "batch_shape", "entropy", "event_shape", "log_prob",
            "support"])

        loc = utils.force_array(numpy.random.uniform(
            -1, 1, self.full_shape).astype(numpy.float32))
        scale = utils.force_array(numpy.exp(numpy.random.uniform(
            -1, 1, self.full_shape)).astype(numpy.float32))

        if self.reinterpreted_batch_ndims is None:
            reinterpreted_batch_ndims = max(0, len(self.inner_shape) - 1)
        else:
            reinterpreted_batch_ndims = self.reinterpreted_batch_ndims

        batch_ndim = len(self.inner_shape) - reinterpreted_batch_ndims
        self.shape = self.inner_shape[:batch_ndim]
        self.event_shape = \
            self.inner_shape[batch_ndim:] + self.inner_event_shape
        d = functools.reduce(operator.mul, self.event_shape, 1)

        if self.event_shape == ():
            self.scipy_dist = stats.norm

            self.params = {"loc": loc, "scale": scale}
            self.scipy_params = {"loc": loc, "scale": scale}

        else:
            self.scipy_dist = stats.multivariate_normal

            scale_tril = numpy.eye(d).astype(numpy.float32) * \
                scale.reshape(self.shape + (d,))[..., None]
            cov = numpy.einsum('...ij,...jk->...ik', scale_tril, scale_tril)

            self.params = {"loc": loc, "scale": scale}
            self.scipy_params = {"mean": numpy.reshape(
                loc, self.shape + (d,)), "cov": cov} 

def check_initializer_statistics(self, backend_config, n):
        from scipy import stats

        xp = backend_config.xp
        ws = numpy.empty((n,) + self.shape, dtype=self.dtype)
        ws = backend_config.get_array(ws)
        for i in range(n):
            initializer = self.target(**self.target_kwargs)
            initializer(xp.squeeze(ws[i:i+1], axis=0))

        fan = self.fan_option or default_fan.get(self.target)
        expected_std = self.scale or default_scale.get(self.target) or 1.
        expected_std *= default_coeff.get(self.target) or 1.
        if fan is not None:
            if fan == 'fan_in':
                expected_std *= math.sqrt(1. / self.fans[0])
            elif fan == 'fan_out':
                expected_std *= math.sqrt(1. / self.fans[1])
            elif fan == 'fan_avg':
                expected_std *= math.sqrt(2. / sum(self.fans))
            else:
                assert False

        sampless = cuda.to_cpu(ws.reshape(n, -1).T)
        alpha = 0.01 / len(sampless)
        for samples in sampless:
            _, p = stats.kstest(samples, stats.norm(0, expected_std).cdf)
            assert p >= alpha 

def ppf(self, u):
        """
        Note: if dim(u) < self.dim, the remaining columns are filled with 0
        Useful in case the distribution is partly degenerate
        """
        N, du = u.shape
        if du < self.dim:
            z = np.zeros((N, self.dim))
            z[:, :du] = stats.norm.ppf(u)
        else:
            z = stats.norm.ppf(u)
        return self.linear_transform(z) 

def rvs(self, size=None):
        if size is None:
            sh = np.broadcast(self.loc, self.scale).shape
            # sh=() when both loc and scale are scalars
            N = 1 if len(sh) == 0 else sh[0]
        else:
            N = size
        z = stats.norm.rvs(size=(N, self.dim))
        return self.linear_transform(z) 

def ppf(self, u):
        return stats.norm.ppf(u, loc=self.loc, scale=self.scale) 

def logpdf(self, x):
        return stats.norm.logpdf(x, loc=self.loc, scale=self.scale) 

def multidmodel():
    """Multidimensional model with normal prior."""
    
    mu = np.array([-6.6, 3, 1.0, -.12])
    sd = np.array([.13, 5, .9, 1.0])
    
    x = SampledParam(norm, loc=mu, scale=sd)
    like = simple_likelihood
    
    return [x], like 

def onedmodel():
    """One dimensional model with normal prior."""
    
    mu = -2
    sd = 3
    x = SampledParam(norm, loc=mu, scale=sd)
    like = simple_likelihood   
    
    return [x], like 

def test_extra_kwarg(self):
        # **kwargs to _pdf are ignored.
        # this is a limitation of the framework (_pdf(x, *goodargs))
        class _distr_gen(stats.rv_continuous):
            def _pdf(self, x, *args, **kwargs):
                # _pdf should handle *args, **kwargs itself.  Here "handling" is
                # ignoring *args and looking for ``extra_kwarg`` and using that.
                extra_kwarg = kwargs.pop('extra_kwarg', 1)
                return stats.norm._pdf(x) * extra_kwarg

        dist = _distr_gen(shapes='extra_kwarg')
        assert_equal(dist.pdf(1, extra_kwarg=3), stats.norm.pdf(1)) 

def test_star_args_2(self):
        # test _pdf with named & starargs
        # NB: **kwargs of pdf will never reach _pdf
        class _dist_gen(stats.rv_continuous):
            def _pdf(self, x, offset, *args):
                extra_kwarg = args[0]
                return stats.norm._pdf(x) * extra_kwarg + offset

        dist = _dist_gen(shapes='offset, extra_kwarg')
        assert_equal(dist.pdf(0.5, offset=111, extra_kwarg=33),
                     stats.norm.pdf(0.5)*33 + 111)
        assert_equal(dist.pdf(0.5, 111, 33),
                     stats.norm.pdf(0.5)*33 + 111)