Standard Deviation
 http://en.wikipedia.org/wiki/Standard_deviation
431337
def standard_deviation(values)
      raise "Values must be an array" unless values.is_a?(Array)

      Math.sqrt(variance(values))
    end
282451
function() {
			// Get the x and y coordinates
			var xCoords = this.pluck('x');
			var yCoords = this.pluck('y');

			// Get the means for the x and y coordinates
			var meanX = xCoords.mean();
			var meanY = yCoords.mean();

			// Get the standard deviations for the x and y coordinates
			var stdDevX = xCoords.stdDev();
			var stdDevY = yCoords.stdDev();

			// Map each element to the difference of the element and the mean
			// divided by the standard deviation
			xCoords.map(function(num) {
				return (num - meanX) / stdDevX;
			});
			yCoords.map(function(num) {
				return (num - meanY) / stdDevY;
			});

			// Multiply each element in the x by the corresponding value in
			// the y
			var nums = this.clone().map(function(num, index) {
				return xCoords.get(index) * yCoords.get(index);
			});

			// r is the sum of xCoords over the number of points minus 1
			return nums.sum() / (nums.size() - 1);
		}
475775
def stddev(self):
    """"""
    if len(self) < 2:
      return float('NaN')
    # The stupidest algorithm, but it works fine.
    try:
      arr = self.samples()
      mean = sum(arr) / len(arr)
      bigsum = 0.0
      for x in arr:
        bigsum += (x - mean)**2
      return sqrt(bigsum / (len(arr) - 1))
    except ZeroDivisionError:
      return float('NaN')
947035
def _synDelayParams(self):
        '''
        
        '''
        delays = {}
        #mean delays
        loc = np.zeros((len(self.y), len(self.X)))
        loc[:, 0] = self.delays[0]
        loc[:, 1::2] = self.delays[0]
        loc[:, 2::2] = self.delays[1]
        #standard deviations
        scale = loc * self.delay_rel_sd
        
        #prepare output
        delay_loc = {}
        for i, y in enumerate(self.y):
            delay_loc.update({y : loc[i]})
        
        delay_scale = {}
        for i, y in enumerate(self.y):
            delay_scale.update({y : scale[i]})
                
        return delay_loc, delay_scale
669550
def _get_stddevs(self, C, stddev_types, num_sites):
        
        assert all(stddev_type in self.DEFINED_FOR_STANDARD_DEVIATION_TYPES
                   for stddev_type in stddev_types)

        if self.DEFINED_FOR_INTENSITY_MEASURE_COMPONENT == 'Random horizontal':
            # Using equation 8 of the paper,
            # corrected as indicated in the erratum
            Sr = np.sqrt(C['SlZ']**2 + (C['S3']/np.sqrt(2))**2)
        else:
            Sr = C['SlZ']

        stddevs = [np.zeros(num_sites) + Sr for _ in stddev_types]
        return stddevs
685535
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        

        sites.vs30 = 700 * np.ones(len(sites.vs30))

        mean, stddevs = super().get_mean_and_stddevs(
            sites, rup, dists, imt, stddev_types)

        C = CauzziFaccioli2008SWISS01.COEFFS

        tau_ss = 'tau'
        log_phi_ss = np.log(10)
        mean, stddevs = _apply_adjustments(
            C, self.COEFFS_FS_ROCK[imt], tau_ss,
            mean, stddevs, sites, rup, dists.rhypo, imt, stddev_types,
            log_phi_ss)

        return mean, np.log(10 ** np.array(stddevs))
1325234
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        assert all(stddev_type in self.DEFINED_FOR_STANDARD_DEVIATION_TYPES
                   for stddev_type in stddev_types)

        C = self.COEFFS[imt]
        mean = self._get_mean(
            C, rup.mag, rup.rake, rup.dip, dists.rrup, dists.rjb
        )
        stddevs = self._get_stddevs(C, rup.mag, stddev_types, dists.rrup.size)

        return mean, stddevs
1085939
def _get_stddevs(self, C, stddev_types, stddev_shape):
        
        stddevs = []
        for stddev_type in stddev_types:
            assert stddev_type in self.DEFINED_FOR_STANDARD_DEVIATION_TYPES
            if stddev_type == const.StdDev.TOTAL:
                stddevs.append(C["sigtot"] + np.zeros(stddev_shape))
            elif stddev_type == const.StdDev.INTRA_EVENT:
                stddevs.append(C['sig2'] + np.zeros(stddev_shape))
            elif stddev_type == const.StdDev.INTER_EVENT:
                stddevs.append(C['sig1'] + np.zeros(stddev_shape))
        return stddevs
1183718
def get_stddevs(self, mean, stddev_types):
        
        mean = np.exp(mean)
        sigma = 0.39 + np.zeros(mean.shape)
        sigma[mean < 0.068] = 0.55
        idx = np.logical_and(mean >= 0.068, mean <= 0.21)
        sigma[idx] = 0.173- 0.140 * np.log(mean[idx])
        stddevs = []
        for stddev in stddev_types:
            if stddev == const.StdDev.TOTAL:
                stddevs.append(sigma)
        return stddevs
378070
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        C = self.COEFFS[imt]
        # Get the mean
        mean = self.get_mean(C, rup, sites, dists)
        if imt.name == "PGV":
            # Convert from log10 m/s to ln cm/s
            mean = np.log((10.0 ** mean) * 100.)
        else:
            # convert from log10 m/s/s to ln g
            mean = np.log((10.0 ** mean) / g)

        # Get the standard deviations
        stddevs = self.get_stddevs(C, mean.shape, stddev_types)
        return mean, stddevs
446669
def _get_stddevs(self, C, stddev_types, rup, imt, num_sites):
        
        stddevs = []
        for stddev_type in stddev_types:
            sigma_mean = self._compute_standard_dev(rup, imt, C)
            sigma_tot = np.sqrt((sigma_mean ** 2) + (C['SigmaReg'] ** 2))
            sigma_tot = np.log10(np.exp(sigma_tot))
            stddevs.append(sigma_tot + np.zeros(num_sites))
        return stddevs
883791
def _get_stddevs(self, sites, rup, C, stddev_types, ln_y_ref, exp1, exp2):
        
        ret = []
        for stddev_type in stddev_types:
            assert stddev_type in self.DEFINED_FOR_STANDARD_DEVIATION_TYPES
            if stddev_type == const.StdDev.TOTAL:
                # eq. 13
                ret.append(0.65 * np.ones_like(sites.vs30))
        return ret
960558
def _get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        # distances
        distsl = copy.copy(dists)
        distsl.rjb, distsl.rrup = \
            utils.get_equivalent_distances_east(rup.mag, dists.repi)
        #
        # Pezeshk et al. 2011 - Rrup
        mean1, stds1 = super().get_mean_and_stddevs(sites, rup, distsl, imt,
                                                    stddev_types)
        mean1 = self.apply_correction_to_BC(mean1, imt, distsl)
        #
        # Atkinson 2008 - Rjb
        gmpe = Atkinson2008prime()
        mean2, stds2 = gmpe.get_mean_and_stddevs(sites, rup, distsl, imt,
                                                 stddev_types)
        #
        # Silva et al. 2002 - Rjb
        gmpe = SilvaEtAl2002SingleCornerSaturation()
        mean4, stds4 = gmpe.get_mean_and_stddevs(sites, rup, distsl, imt,
                                                 stddev_types)
        mean4 = self.apply_correction_to_BC(mean4, imt, distsl)
        #
        # Silva et al. 2002 - Rjb
        gmpe = SilvaEtAl2002DoubleCornerSaturation()
        mean5, stds5 = gmpe.get_mean_and_stddevs(sites, rup, distsl, imt,
                                                 stddev_types)
        mean5 = self.apply_correction_to_BC(mean5, imt, distsl)
        #
        # distances
        distsl.rjb, distsl.rrup = \
            utils.get_equivalent_distances_east(rup.mag, dists.repi, ab06=True)
        #
        # Atkinson and Boore 2006 - Rrup
        gmpe = AtkinsonBoore2006Modified2011()
        mean3, stds3 = gmpe.get_mean_and_stddevs(sites, rup, distsl, imt,
                                                 stddev_types)
        # Computing adjusted mean and stds
        mean_adj = mean1*0.2 + mean2*0.2 + mean3*0.2 + mean4*0.2 + mean5*0.2

        # Note that in this case we do not apply a triangular smoothing on
        # distance as explained at page 996 of Atkinson and Adams (2013)
        # for the calculation of the standard deviation
        stds_adj = np.log(np.exp(stds1)*0.2 + np.exp(stds2)*0.2 +
                          np.exp(stds3)*0.2 + np.exp(stds4)*0.2 +
                          np.exp(stds5)*0.2)
        #
        return mean_adj, stds_adj
473455
def _get_stddevs(self, C, mag, stddev_types, sites):
        
        num_sites = sites.vs30.size
        sigma_intra = np.zeros(num_sites)

        # interevent stddev
        tau = sigma_intra + C['tau']

        # intraevent std (equations 8a-8c page 29)
        if mag < 5.0:
            sigma_intra += C['sigmaM6'] - C['sigSlope']
        elif 5.0 <= mag < 7.0:
            sigma_intra += C['sigmaM6'] + C['sigSlope'] * (mag - 6)
        else:
            sigma_intra += C['sigmaM6'] + C['sigSlope']

        std = []

        for stddev_type in stddev_types:
            if stddev_type == const.StdDev.TOTAL:
                # equation 9 page 29
                std += [np.sqrt(sigma_intra**2 + tau**2)]
            elif stddev_type == const.StdDev.INTRA_EVENT:
                std.append(sigma_intra)
            elif stddev_type == const.StdDev.INTER_EVENT:
                std.append(tau)

        return std
1312455
def _sd_of_runs(stats, mean, key='runs'):
    

    num_runs = len(stats[key])
    first = stats[key][0]

    standard_deviation = {}
    for stat_key in first:
        # Skip non numberic attributes
        if isinstance(first[stat_key], numbers.Number):
            standard_deviation[stat_key] = math.sqrt(
                sum((run[stat_key] - mean[stat_key])**2
                    for run in stats[key]) / float(num_runs))

    return standard_deviation
1030233
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        C = self.COEFFS[imt]

        mean = (self._compute_magnitude_term(C, rup.mag) +
                self._compute_distance_term(C, dists.rhypo, rup.mag))
        stddevs = self._get_stddevs(C, dists.rhypo, stddev_types)
        return mean, stddevs
1048914
def stats():
    ''''''
    import re
    import sys
    import math
    values = []
    for line in sys.stdin:
        values.extend(map(float, re.findall(r'\d+\.?\d+', line)))

    mean = sum(values) / len(values)
    variance = sum((val - mean) ** 2 for val in values) / len(values)
    print '%3i items; mean: %10.5f; std-dev: %10.5f' % (
        len(values), mean, math.sqrt(variance))