我们从Python开源项目中，提取了以下22个代码示例，用于说明如何使用scipy.special.binom()。
def embed(self, G, rng): if G.order() < self.cycle_length: return np.array([0.0], dtype = "f") succ = 0 max_cycle = sp.binom(G.order(), self.cycle_length) for _ in range(self.sample_size): us = rng.choice(G.nodes(), self.cycle_length) H = G.subgraph(us.tolist()) for c in nx.simple_cycles(H): if len(c) == self.cycle_length: succ += 1 break val = [max_cycle * (succ / self.sample_size)] return np.array(val, dtype = "f")
def embed(self, G, rng): if G.order() < self.cycle_length: return np.array([0.0], dtype = "f") succ, samples = 0, 0 max_cycle = sp.binom(G.order(), self.cycle_length) for _ in range(self.sample_cap): us = rng.choice(G.nodes(), self.cycle_length) H = G.subgraph(us.tolist()) samples += 1 for c in nx.simple_cycles(H): if len(c) == self.cycle_length: succ += 1 break if succ >= self.successes: break return np.array([max_cycle * (succ / samples)], dtype = "f")
def test_beta_binomial_two_identical_models(db_path, sampler): binomial_n = 5 def model_fun(args): return {"result": st.binom(binomial_n, args.theta).rvs()} models = [model_fun for _ in range(2)] models = list(map(SimpleModel, models)) population_size = ConstantPopulationSize(800) parameter_given_model_prior_distribution = [Distribution(theta=st.beta( 1, 1)) for _ in range(2)] abc = ABCSMC(models, parameter_given_model_prior_distribution, MinMaxDistanceFunction(measures_to_use=["result"]), population_size, eps=MedianEpsilon(.1), sampler=sampler) abc.new(db_path, {"result": 2}) minimum_epsilon = .2 history = abc.run(minimum_epsilon, max_nr_populations=3) mp = history.get_model_probabilities(history.max_t) assert abs(mp.p[0] - .5) + abs(mp.p[1] - .5) < .08
def test_beta_binomial_two_identical_models_adaptive(db_path, sampler): binomial_n = 5 def model_fun(args): return {"result": st.binom(binomial_n, args.theta).rvs()} models = [model_fun for _ in range(2)] models = list(map(SimpleModel, models)) population_size = AdaptivePopulationSize(800) parameter_given_model_prior_distribution = [ Distribution(theta=st.beta(1, 1)) for _ in range(2)] abc = ABCSMC(models, parameter_given_model_prior_distribution, MinMaxDistanceFunction(measures_to_use=["result"]), population_size, eps=MedianEpsilon(.1), sampler=sampler) abc.new(db_path, {"result": 2}) minimum_epsilon = .2 history = abc.run(minimum_epsilon, max_nr_populations=3) mp = history.get_model_probabilities(history.max_t) assert abs(mp.p[0] - .5) + abs(mp.p[1] - .5) < .08
def shift_polynomial(coef, d): # given polynomial sum_k a_k x^k -> sum_k a_k (x+d)^k coef2 = np.zeros(5, dtype="float64") for i in range(5): for j in range(i, 5): coef2[i] = coef2[i] + coef[j] * (d**(j - i)) * sp.binom(j, i) return coef2 # def shift_graph_2(graph1, graph2, d): # # numpy matrix implementation of shift_graph # A = np.zeros((5,5),dtype="float64") # for i in range(5): # for j in range(i,5): # A[j,i] = (d**(j-i)) * sp.binom(j,i) # graph2[:,3:] = graph2[:,3:].dot(A) # return graph2
def expected_fscore(n_predicted, n=100, prob_per_instance=.01): #assert n_predicted in (1, n) ppi = prob_per_instance ppibar = (1 - prob_per_instance) def fscore(fn, tp, fp): num = 2*tp denom = 2*tp + fp + fn return num / denom parts = [] for ntrue in range(1, n+1): # ntrue := this many were actually true in the gold standard if n_predicted == 1: tp = 1 fn = ntrue - 1 fp = 0 # Prob. of our 1 predicted label being true prob_weight = ppi # Prob of ntrue-1 other labels being true out of n-1 prob_weight *= binom(n-1, ntrue-1) * ppi**(ntrue-1) * ppibar**(n-ntrue) elif n_predicted == n: tp = ntrue fn = 0 fp = n - ntrue prob_weight = binom(n, ntrue) * ppi**ntrue * ppibar**(n-ntrue) else: assert False, "fuck it" for tp in range(1, ntrue): fp = n_predicted - tp fn = ntrue - tp fs = fscore(fn, tp, fp) dbg('fscore = {}'.format(fs)) dbg('prob weight = {}'.format(prob_weight)) res = fs * prob_weight parts.append(res) return parts
def calculate_initial_centroids_first_run( baskets, neighbors, no_neighbors, item_baskets, use_neighbors, num_rnd_init, verbose=False): max_no_nbr = list() max_no_nbr_index = 0 if use_neighbors: for b in baskets: no_nbr_list = no_neighbors[b] no_nbr_tuple = (b, no_nbr_list) no_nbr = no_nbr_list.count() if no_nbr > max_no_nbr_index: max_no_nbr_index = no_nbr max_no_nbr = list() max_no_nbr.append(no_nbr_tuple) elif no_nbr == max_no_nbr_index: max_no_nbr.append(no_nbr_tuple) if max_no_nbr_index == 0 or len(max_no_nbr) == len(baskets): # print 'rnd' if verbose: print datetime.datetime.now(), 'init random' num_rnd_init = min(num_rnd_init, int(binom(len(baskets), 2))) m0_index, m1_index = random_init(baskets, neighbors, use_neighbors, num=num_rnd_init) else: # print 'no nbr' if verbose: print datetime.datetime.now(), 'init no neighbors' m0_index, m1_index = no_nbr_init(max_no_nbr, baskets, neighbors, use_neighbors) return m0_index, m1_index
def size(N, nb): """ Return the size of the boson many-body basis. Parameters ---------- N : int Number of sites nb : int Number of bosons """ return int(binom(nb + N - 1, nb))
def gen_probs(n, a, b): probs = np.zeros(n+1) for k in range(n+1): probs[k] = binom(n, k) * beta(k + a, n - k + b) / beta(a, b) return probs
def genton(X): """ Return the Genton Variogram of the given sample X. X has to be an even-length array of point pairs like: x1, x1+h, x2, x2+h ...., xn, xn + h. If X.ndim > 1, genton will be called recursively and a list of Cressie-Hawkins Variances is returned. Genton, M. G., (1998): Highly robust variogram estimation, Math. Geol., 30, 213 - 221. :param X: :return: """ _X = np.array(X) if any([isinstance(_, list) or isinstance(_, np.ndarray) for _ in _X]): return np.array([genton(_) for _ in _X]) # check even if len(_X) % 2 > 0: raise ValueError('The sample does not have an even length: {}'.format(_X)) # calculate try: y = [np.abs( (_X[i] - _X[i + 1]) - (_X[j]) - _X[j + 1] ) for i in np.arange(0, len(_X), 2) for j in np.arange(0, len(_X), 2) if i < j] # get k k is binom(N(x)/2+1, 2) k = binom(int(len(_X) / 2 + 1), 2) # return the kth percentile return 0.5 * np.power(2.219 * np.percentile(y, k), 2) except ZeroDivisionError: return np.nan
def binomial(n): ''' Return all binomial coefficents for a given order. For n > 5, scipy.special.binom is used, below we hardcode to avoid the scipy.special dependancy. ''' if n == 1: return [1,1] elif n == 2: return [1,2,1] elif n == 3: return [1,3,3,1] elif n == 4: return [1,4,6,4,1] elif n == 5: return [1,5,10,10,5,1] else: from scipy.special import binom return binom(n,np.arange(n+1))
def _F_mk(Q_p, P, n_p, J, r_b, r_p, z, pikg, sigma): """ Complex matrix F_mk from Claesson and Hellstrom (2011), EQ. 34 """ F = np.zeros((n_p, J), dtype=np.cfloat) for m in range(n_p): for k in range(J): fmk = 0. + 0.j for n in range(n_p): # First term if m != n: fmk += Q_p[n]*pikg/(k+1)*(r_p[m]/(z[n] - z[m]))**(k+1) # Second term fmk += sigma*Q_p[n]*pikg/(k+1)*(r_p[m]*np.conj(z[n])/( r_b**2 - z[m]*np.conj(z[n])))**(k+1) for j in range(J): # Third term if m != n: fmk += P[n,j]*binom(j+k+1, j) \ *r_p[n]**(j+1)*(-r_p[m])**(k+1) \ /(z[m] - z[n])**(j+k+2) # Fourth term j_pend = np.min((k, j)) + 1 for jp in range(j_pend): fmk += sigma*np.conj(P[n,j])*binom(j+1, jp) \ *binom(j+k-jp+1, j)*r_p[n]**(j+1) \ *r_p[m]**(k+1)*z[m]**(j+1-jp) \ *np.conj(z[n])**(k+1-jp) \ /(r_b**2 - z[m]*np.conj(z[n]))**(k+j+2-jp) F[m,k] = fmk return F
def readout(q, r): """C matrix to decode a delay of r*theta from the delay state for theta. ``r`` is a ratio between 0 (``t=0``) and 1 (``t=-theta``). """ c = np.zeros(q) for i in range(q): j = np.arange(i+1, dtype=np.float64) c[q-1-i] += 1 / binom(q, i) * np.sum( binom(q, j) * binom(2*q - 1 - j, i - j) * (-r)**(i - j)) return c
def test_beta_binomial_different_priors(db_path, sampler): binomial_n = 5 def model(args): return {"result": st.binom(binomial_n, args['theta']).rvs()} models = [model for _ in range(2)] models = list(map(SimpleModel, models)) population_size = ConstantPopulationSize(800) a1, b1 = 1, 1 a2, b2 = 10, 1 parameter_given_model_prior_distribution = [Distribution(theta=RV("beta", a1, b1)), Distribution(theta=RV("beta", a2, b2))] abc = ABCSMC(models, parameter_given_model_prior_distribution, MinMaxDistanceFunction(measures_to_use=["result"]), population_size, eps=MedianEpsilon(.1), sampler=sampler) n1 = 2 abc.new(db_path, {"result": n1}) minimum_epsilon = .2 history = abc.run(minimum_epsilon, max_nr_populations=3) mp = history.get_model_probabilities(history.max_t) def B(a, b): return gamma(a) * gamma(b) / gamma(a + b) def expected_p(a, b, n1): return binom(binomial_n, n1) * B(a + n1, b + binomial_n - n1) / B(a, b) p1_expected_unnormalized = expected_p(a1, b1, n1) p2_expected_unnormalized = expected_p(a2, b2, n1) p1_expected = p1_expected_unnormalized / (p1_expected_unnormalized + p2_expected_unnormalized) p2_expected = p2_expected_unnormalized / (p1_expected_unnormalized + p2_expected_unnormalized) assert abs(mp.p[0] - p1_expected) + abs(mp.p[1] - p2_expected) < .08
def shift_polynomial(coef, d): # given polynomial sum_k a_k x^k -> sum_k a_k (x+d)^k coef2 = np.zeros(5, dtype="float64") for i in range(5): for j in range(i, 5): coef2[i] = coef2[i] + coef[j] * (d**(j - i)) * sp.binom(j, i) return coef2
def _check_subparams(self, n_samples, n_features): n_subsamples = self.n_subsamples if self.fit_intercept: n_dim = n_features + 1 else: n_dim = n_features if n_subsamples is not None: if n_subsamples > n_samples: raise ValueError("Invalid parameter since n_subsamples > " "n_samples ({0} > {1}).".format(n_subsamples, n_samples)) if n_samples >= n_features: if n_dim > n_subsamples: plus_1 = "+1" if self.fit_intercept else "" raise ValueError("Invalid parameter since n_features{0} " "> n_subsamples ({1} > {2})." "".format(plus_1, n_dim, n_samples)) else: # if n_samples < n_features if n_subsamples != n_samples: raise ValueError("Invalid parameter since n_subsamples != " "n_samples ({0} != {1}) while n_samples " "< n_features.".format(n_subsamples, n_samples)) else: n_subsamples = min(n_dim, n_samples) if self.max_subpopulation <= 0: raise ValueError("Subpopulation must be strictly positive " "({0} <= 0).".format(self.max_subpopulation)) all_combinations = max(1, np.rint(binom(n_samples, n_subsamples))) n_subpopulation = int(min(self.max_subpopulation, all_combinations)) return n_subsamples, n_subpopulation
def point_probability(k, n, l, p=0.5): return binom(n - l, k) * p ** k * (1 - p) ** (n - k - l)
def testLogCombinationsBinomial(self): n = [2, 5, 12, 15] k = [1, 2, 4, 11] log_combs = np.log(special.binom(n, k)) with self.test_session(): n = np.array(n, dtype=np.float32) counts = [[1., 1], [2., 3], [4., 8], [11, 4]] log_binom = distribution_util.log_combinations(n, counts) self.assertEqual([4], log_binom.get_shape()) self.assertAllClose(log_combs, log_binom.eval())
def sample(self, N=100): """Compute N points of the Bezier curve. Returns two lists, containing the x and y coordinates.""" def p(t, coord): n = len(coord) res = 0 for k in range(n): res += coord[k] * binom(n - 1, k) * t ** k * (1 - t) ** (n - 1 - k) return (res) x = [p(t / N, self.controlPoints_x) for t in range(N)] y = [p(t / N, self.controlPoints_y) for t in range(N)] return x, y
def estimation(self, y, ds=None): """ Estimate the fourth multivariate extension of Spearman's rho. Parameters ---------- y : (number of samples, dimension)-ndarray One row of y corresponds to one sample. ds : int vector, vector of ones ds[i] = 1 (for all i): the i^th subspace is one-dimensional. If ds is not given (ds=None), the vector of ones [ds = ones(y.shape[1],dtype='int')] is emulated inside the function. Returns ------- a : float Estimated fourth multivariate extension of Spearman's rho. References ---------- Friedrich Shmid, Rafael Schmidt, Thomas Blumentritt, Sandra Gaiser, and Martin Ruppert. Copula Theory and Its Applications, Chapter Copula based Measures of Multivariate Association. Lecture Notes in Statistics. Springer, 2010. Friedrich Schmid and Rafael Schmidt. Multivariate extensions of Spearman's rho and related statistics. Statistics & Probability Letters, 77:407-416, 2007. Maurice G. Kendall. Rank correlation methods. London, Griffin, 1970. C. Spearman. The proof and measurement of association between two things. The American Journal of Psychology, 15:72-101, 1904. Examples -------- a = co.estimation(y,ds) """ if ds is None: # emulate 'ds = vector of ones' ds = ones(y.shape[1], dtype='int') # verification: self.verification_compatible_subspace_dimensions(y, ds) self.verification_one_dimensional_subspaces(ds) num_of_samples, dim = y.shape # number of samples, dimension u = copula_transformation(y) m_triu = triu(ones((dim, dim)), 1) # upper triangular mask b = binom(dim, 2) a = 12 * sum(dot((1 - u).T, (1 - u)) * m_triu) /\ (b * num_of_samples) - 3 return a
def test_beta_binomial_different_priors_initial_epsilon_from_sample(db_path, sampler): binomial_n = 5 def model(args): return {"result": st.binom(binomial_n, args.theta).rvs()} models = [model for _ in range(2)] models = list(map(SimpleModel, models)) population_size = ConstantPopulationSize(800) a1, b1 = 1, 1 a2, b2 = 10, 1 parameter_given_model_prior_distribution = [Distribution(theta=RV("beta", a1, b1)), Distribution(theta=RV("beta", a2, b2))] abc = ABCSMC(models, parameter_given_model_prior_distribution, MinMaxDistanceFunction(measures_to_use=["result"]), population_size, eps=MedianEpsilon(median_multiplier=.9), sampler=sampler) n1 = 2 abc.new(db_path, {"result": n1}) minimum_epsilon = -1 history = abc.run(minimum_epsilon, max_nr_populations=5) mp = history.get_model_probabilities(history.max_t) def B(a, b): return gamma(a) * gamma(b) / gamma(a + b) def expected_p(a, b, n1): return binom(binomial_n, n1) * B(a + n1, b + binomial_n - n1) / B(a, b) p1_expected_unnormalized = expected_p(a1, b1, n1) p2_expected_unnormalized = expected_p(a2, b2, n1) p1_expected = p1_expected_unnormalized / (p1_expected_unnormalized + p2_expected_unnormalized) p2_expected = p2_expected_unnormalized / (p1_expected_unnormalized + p2_expected_unnormalized) assert abs(mp.p[0] - p1_expected) + abs(mp.p[1] - p2_expected) < .08
def fit(self, X, y): """Fit linear model. Parameters ---------- X : numpy array of shape [n_samples, n_features] Training data y : numpy array of shape [n_samples] Target values Returns ------- self : returns an instance of self. """ random_state = check_random_state(self.random_state) X, y = check_X_y(X, y, y_numeric=True) n_samples, n_features = X.shape n_subsamples, self.n_subpopulation_ = self._check_subparams(n_samples, n_features) self.breakdown_ = _breakdown_point(n_samples, n_subsamples) if self.verbose: print("Breakdown point: {0}".format(self.breakdown_)) print("Number of samples: {0}".format(n_samples)) tol_outliers = int(self.breakdown_ * n_samples) print("Tolerable outliers: {0}".format(tol_outliers)) print("Number of subpopulations: {0}".format( self.n_subpopulation_)) # Determine indices of subpopulation if np.rint(binom(n_samples, n_subsamples)) <= self.max_subpopulation: indices = list(combinations(range(n_samples), n_subsamples)) else: indices = [choice(n_samples, size=n_subsamples, replace=False, random_state=random_state) for _ in range(self.n_subpopulation_)] n_jobs = _get_n_jobs(self.n_jobs) index_list = np.array_split(indices, n_jobs) weights = Parallel(n_jobs=n_jobs, verbose=self.verbose)( delayed(_lstsq)(X, y, index_list[job], self.fit_intercept) for job in range(n_jobs)) weights = np.vstack(weights) self.n_iter_, coefs = _spatial_median(weights, max_iter=self.max_iter, tol=self.tol) if self.fit_intercept: self.intercept_ = coefs[0] self.coef_ = coefs[1:] else: self.intercept_ = 0. self.coef_ = coefs return self
