我们从Python开源项目中，提取了以下50个代码示例，用于说明如何使用numpy.random.rand()。
def test_basic(self): y1 = np.array([1, 2, 3]) assert_(average(y1, axis=0) == 2.) y2 = np.array([1., 2., 3.]) assert_(average(y2, axis=0) == 2.) y3 = [0., 0., 0.] assert_(average(y3, axis=0) == 0.) y4 = np.ones((4, 4)) y4[0, 1] = 0 y4[1, 0] = 2 assert_almost_equal(y4.mean(0), average(y4, 0)) assert_almost_equal(y4.mean(1), average(y4, 1)) y5 = rand(5, 5) assert_almost_equal(y5.mean(0), average(y5, 0)) assert_almost_equal(y5.mean(1), average(y5, 1)) y6 = np.matrix(rand(5, 5)) assert_array_equal(y6.mean(0), average(y6, 0))
def test_basic(self): from numpy.random import rand a = rand(20, 10, 10, 1, 1) b = rand(20, 1, 10, 1, 20) c = rand(1, 1, 20, 10) assert_array_equal(np.squeeze(a), np.reshape(a, (20, 10, 10))) assert_array_equal(np.squeeze(b), np.reshape(b, (20, 10, 20))) assert_array_equal(np.squeeze(c), np.reshape(c, (20, 10))) # Squeezing to 0-dim should still give an ndarray a = [[[1.5]]] res = np.squeeze(a) assert_equal(res, 1.5) assert_equal(res.ndim, 0) assert_equal(type(res), np.ndarray)
def mi(x, y, k=3, base=2): """ Mutual information of x and y x, y should be a list of vectors, e.g. x = [[1.3], [3.7], [5.1], [2.4]] if x is a one-dimensional scalar and we have four samples """ assert len(x) == len(y), "Lists should have same length" assert k <= len(x) - 1, "Set k smaller than num. samples - 1" intens = 1e-10 # small noise to break degeneracy, see doc. x = [list(p + intens * nr.rand(len(x[0]))) for p in x] y = [list(p + intens * nr.rand(len(y[0]))) for p in y] points = zip2(x, y) # Find nearest neighbors in joint space, p=inf means max-norm tree = ss.cKDTree(points) dvec = [tree.query(point, k + 1, p=float('inf'))[0][k] for point in points] a, b, c, d = avgdigamma(x, dvec), avgdigamma(y, dvec), digamma(k), digamma(len(x)) return (-a - b + c + d) / log(base)
def cmi(x, y, z, k=3, base=2): """ Mutual information of x and y, conditioned on z x, y, z should be a list of vectors, e.g. x = [[1.3], [3.7], [5.1], [2.4]] if x is a one-dimensional scalar and we have four samples """ assert len(x) == len(y), "Lists should have same length" assert k <= len(x) - 1, "Set k smaller than num. samples - 1" intens = 1e-10 # small noise to break degeneracy, see doc. x = [list(p + intens * nr.rand(len(x[0]))) for p in x] y = [list(p + intens * nr.rand(len(y[0]))) for p in y] z = [list(p + intens * nr.rand(len(z[0]))) for p in z] points = zip2(x, y, z) # Find nearest neighbors in joint space, p=inf means max-norm tree = ss.cKDTree(points) dvec = [tree.query(point, k + 1, p=float('inf'))[0][k] for point in points] a, b, c, d = avgdigamma(zip2(x, z), dvec), avgdigamma(zip2(y, z), dvec), avgdigamma(z, dvec), digamma(k) return (-a - b + c + d) / log(base)
def gen_graphing_report_fn (pct, prefix="gibbs-dist"): def report_fn (net, arrv, gi): global i_ggrf for evt in arrv: if random.rand() < pct: # graph it f = open ("%s%d.txt" % (prefix, i_ggrf), "w") dfn = generate_true_dfn (net, arrv, evt) L,U = dfn_range (net, arrv, evt) eps = (U-L)/100 x = L for i in range(100): f.write ("%.5f %.5f\n" % (x, dfn(x))) x += eps f.close() i_ggrf += 1 return report_fn # For model selection
def sim_likelihood(p1, p2, p3): """Simulates the likelihood.""" # service times (uniformly distributed) sts = (p2 - p1) * rng.rand(n_sim_steps) + p1 # interarrival times (exponentially distributed) iats = -np.log(1.0 - rng.rand(n_sim_steps)) / p3 # arrival times ats = np.cumsum(iats) # interdeparture and departure times idts = np.empty(n_sim_steps) dts = np.empty(n_sim_steps) idts[0] = sts[0] + ats[0] dts[0] = idts[0] for i in xrange(1, n_sim_steps): idts[i] = sts[i] + max(0.0, ats[i] - dts[i-1]) dts[i] = dts[i-1] + idts[i] return sts, iats, ats, idts, dts
def discrete_sample(p, n_samples=1): """ Samples from a discrete distribution. :param p: a distribution with N elements :param n_samples: number of samples :return: vector of samples """ # check distribution #assert isdistribution(p), 'Probabilities must be non-negative and sum to one.' # cumulative distribution c = np.cumsum(p[:-1])[np.newaxis, :] # get the samples r = rng.rand(n_samples, 1) return np.sum((r > c).astype(int), axis=1)
def test_long_sequence(): N_iters = 200 vect_length = 10 Ns = [] Ms = [] for i in range(N_iters): Ms.append(npr.randint(200) + 1) Ns.append(npr.randint(Ms[-1], size=vect_length)) store = BitStore(vect_length) coinflips = npr.rand(N_iters) new_Ns = [] for N, M, r in zip(Ns, Ms, coinflips): if r < 0.75: store.push(N, M) else: new_Ns.append(store.pop(M)) for N, M, r in zip(Ns, Ms, coinflips)[::-1]: if r < 0.75: cur_N = store.pop(M) assert np.all(cur_N == N) else: store.push(new_Ns.pop(), M)
def createSampleData(m = 100, n = 20, scale = 2, p = 0.5): print("Creating sample data.") data = zeros((m, n)) row = scale * random.rand(n) k = 0 for i in range(m): u = random.rand() if u > p: row = scale * random.rand(n) k += 1 data[i] = row random.shuffle(data) for i in range(m): for j in range(n): data[i, j] = random.poisson(data[i, j]) print("Sample data created with {} different cell types.".format(k)) return data
def mi(x, y, k=3, base=2): """Mutual information of x and y. x,y should be a list of vectors, e.g. x = [[1.3], [3.7], [5.1], [2.4]] if x is a one-dimensional scalar and we have four samples. """ assert len(x)==len(y), 'Lists should have same length.' assert k <= len(x) - 1, 'Set k smaller than num samples - 1.' intens = 1e-10 # Small noise to break degeneracy, see doc. x = [list(p + intens*nr.rand(len(x[0]))) for p in x] y = [list(p + intens*nr.rand(len(y[0]))) for p in y] points = zip2(x,y) # Find nearest neighbors in joint space, p=inf means max-norm. tree = ss.cKDTree(points) dvec = [tree.query(point, k+1, p=float('inf'))[0][k] for point in points] a = avgdigamma(x,dvec) b = avgdigamma(y,dvec) c = digamma(k) d = digamma(len(x)) return (-a-b+c+d) / log(base)
def cmi(x, y, z, k=3, base=2): """Mutual information of x and y, conditioned on z x,y,z should be a list of vectors, e.g. x = [[1.3], [3.7], [5.1], [2.4]] if x is a one-dimensional scalar and we have four samples """ assert len(x)==len(y), 'Lists should have same length.' assert k <= len(x) - 1, 'Set k smaller than num samples - 1.' intens = 1e-10 # Small noise to break degeneracy, see doc. x = [list(p + intens*nr.rand(len(x[0]))) for p in x] y = [list(p + intens*nr.rand(len(y[0]))) for p in y] z = [list(p + intens*nr.rand(len(z[0]))) for p in z] points = zip2(x,y,z) # Find nearest neighbors in joint space, p=inf means max-norm. tree = ss.cKDTree(points) dvec = [tree.query(point, k+1, p=float('inf'))[0][k] for point in points] a = avgdigamma(zip2(x,z), dvec) b = avgdigamma(zip2(y,z), dvec) c = avgdigamma(z,dvec) d = digamma(k) return (-a-b+c+d) / log(base)
def test_np(): npr.seed(0) nx, nineq, neq = 4, 6, 7 Q = npr.randn(nx, nx) G = npr.randn(nineq, nx) A = npr.randn(neq, nx) D = np.diag(npr.rand(nineq)) K_ = np.bmat(( (Q, np.zeros((nx, nineq)), G.T, A.T), (np.zeros((nineq, nx)), D, np.eye(nineq), np.zeros((nineq, neq))), (G, np.eye(nineq), np.zeros((nineq, nineq + neq))), (A, np.zeros((neq, nineq + nineq + neq))) )) K = block(( (Q, 0, G.T, A.T), (0, D, 'I', 0), (G, 'I', 0, 0), (A, 0, 0, 0) )) assert np.allclose(K_, K)
def test_frame_negate(self): expr = self.ex('-') # float lhs = DataFrame(randn(5, 2)) expect = -lhs result = pd.eval(expr, engine=self.engine, parser=self.parser) assert_frame_equal(expect, result) # int lhs = DataFrame(randint(5, size=(5, 2))) expect = -lhs result = pd.eval(expr, engine=self.engine, parser=self.parser) assert_frame_equal(expect, result) # bool doesn't work with numexpr but works elsewhere lhs = DataFrame(rand(5, 2) > 0.5) if self.engine == 'numexpr': with tm.assertRaises(NotImplementedError): result = pd.eval(expr, engine=self.engine, parser=self.parser) else: expect = -lhs result = pd.eval(expr, engine=self.engine, parser=self.parser) assert_frame_equal(expect, result)
def test_series_negate(self): expr = self.ex('-') # float lhs = Series(randn(5)) expect = -lhs result = pd.eval(expr, engine=self.engine, parser=self.parser) assert_series_equal(expect, result) # int lhs = Series(randint(5, size=5)) expect = -lhs result = pd.eval(expr, engine=self.engine, parser=self.parser) assert_series_equal(expect, result) # bool doesn't work with numexpr but works elsewhere lhs = Series(rand(5) > 0.5) if self.engine == 'numexpr': with tm.assertRaises(NotImplementedError): result = pd.eval(expr, engine=self.engine, parser=self.parser) else: expect = -lhs result = pd.eval(expr, engine=self.engine, parser=self.parser) assert_series_equal(expect, result)
def test_arithmetic_interaction(self): index = self.frame.index obj_index = index.asobject dseries = Series(rand(len(index)), index=index) oseries = Series(dseries.values, index=obj_index) result = dseries + oseries expected = dseries * 2 tm.assertIsInstance(result.index, DatetimeIndex) assert_series_equal(result, expected) result = dseries + oseries[:5] expected = dseries + dseries[:5] tm.assertIsInstance(result.index, DatetimeIndex) assert_series_equal(result, expected)
def test_line_lim(self): df = DataFrame(rand(6, 3), columns=['x', 'y', 'z']) ax = df.plot() xmin, xmax = ax.get_xlim() lines = ax.get_lines() self.assertEqual(xmin, lines[0].get_data()[0][0]) self.assertEqual(xmax, lines[0].get_data()[0][-1]) ax = df.plot(secondary_y=True) xmin, xmax = ax.get_xlim() lines = ax.get_lines() self.assertEqual(xmin, lines[0].get_data()[0][0]) self.assertEqual(xmax, lines[0].get_data()[0][-1]) axes = df.plot(secondary_y=True, subplots=True) self._check_axes_shape(axes, axes_num=3, layout=(3, 1)) for ax in axes: self.assertTrue(hasattr(ax, 'left_ax')) self.assertFalse(hasattr(ax, 'right_ax')) xmin, xmax = ax.get_xlim() lines = ax.get_lines() self.assertEqual(xmin, lines[0].get_data()[0][0]) self.assertEqual(xmax, lines[0].get_data()[0][-1])
def test_area_lim(self): df = DataFrame(rand(6, 4), columns=['x', 'y', 'z', 'four']) neg_df = -df for stacked in [True, False]: ax = _check_plot_works(df.plot.area, stacked=stacked) xmin, xmax = ax.get_xlim() ymin, ymax = ax.get_ylim() lines = ax.get_lines() self.assertEqual(xmin, lines[0].get_data()[0][0]) self.assertEqual(xmax, lines[0].get_data()[0][-1]) self.assertEqual(ymin, 0) ax = _check_plot_works(neg_df.plot.area, stacked=stacked) ymin, ymax = ax.get_ylim() self.assertEqual(ymax, 0)
def test_kde_colors(self): tm._skip_if_no_scipy() _skip_if_no_scipy_gaussian_kde() from matplotlib import cm custom_colors = 'rgcby' df = DataFrame(rand(5, 5)) ax = df.plot.kde(color=custom_colors) self._check_colors(ax.get_lines(), linecolors=custom_colors) tm.close() ax = df.plot.kde(colormap='jet') rgba_colors = lmap(cm.jet, np.linspace(0, 1, len(df))) self._check_colors(ax.get_lines(), linecolors=rgba_colors) tm.close() ax = df.plot.kde(colormap=cm.jet) rgba_colors = lmap(cm.jet, np.linspace(0, 1, len(df))) self._check_colors(ax.get_lines(), linecolors=rgba_colors)
def test_partially_invalid_plot_data(self): with tm.RNGContext(42): df = DataFrame(randn(10, 2), dtype=object) df[np.random.rand(df.shape[0]) > 0.5] = 'a' for kind in plotting._common_kinds: if not _ok_for_gaussian_kde(kind): continue with tm.assertRaises(TypeError): df.plot(kind=kind) with tm.RNGContext(42): # area plot doesn't support positive/negative mixed data kinds = ['area'] df = DataFrame(rand(10, 2), dtype=object) df[np.random.rand(df.shape[0]) > 0.5] = 'a' for kind in kinds: with tm.assertRaises(TypeError): df.plot(kind=kind)
def test_errorbar_asymmetrical(self): np.random.seed(0) err = np.random.rand(3, 2, 5) data = np.random.randn(5, 3) df = DataFrame(data) ax = df.plot(yerr=err, xerr=err / 2) self.assertEqual(ax.lines[7].get_ydata()[0], data[0, 1] - err[1, 0, 0]) self.assertEqual(ax.lines[8].get_ydata()[0], data[0, 1] + err[1, 1, 0]) self.assertEqual(ax.lines[5].get_xdata()[0], -err[1, 0, 0] / 2) self.assertEqual(ax.lines[6].get_xdata()[0], err[1, 1, 0] / 2) with tm.assertRaises(ValueError): df.plot(yerr=err.T) tm.close()
def get_grads(nBatch=1, nz=10, neq=1, nineq=3, Qscale=1., Gscale=1., hscale=1., Ascale=1., bscale=1.): assert(nBatch == 1) npr.seed(1) L = np.random.randn(nz, nz) Q = Qscale * L.dot(L.T) G = Gscale * npr.randn(nineq, nz) # h = hscale*npr.randn(nineq) z0 = npr.randn(nz) s0 = npr.rand(nineq) h = G.dot(z0) + s0 A = Ascale * npr.randn(neq, nz) # b = bscale*npr.randn(neq) b = A.dot(z0) p = npr.randn(nBatch, nz) # print(np.linalg.norm(p)) truez = npr.randn(nBatch, nz) Q, p, G, h, A, b, truez = [x.astype(np.float64) for x in [Q, p, G, h, A, b, truez]] _, zhat, nu, lam, slacks = qp_cvxpy.forward_single_np(Q, p[0], G, h, A, b) grads = get_grads_torch(Q, p, G, h, A, b, truez) return [p[0], Q, G, h, A, b, truez], grads
def test_pairwise_jsd_equal_to_jsd_for_pairs_of_perspectives(): nTopics = 4 OT = 7 cn = [str(t) for t in range(nTopics)] perspectives = ['p0', 'p1', 'p2'] # generate random opinions opinions = {} for p in perspectives: o = rand(OT, nTopics) opinions[p] = DataFrame(o / sum(o, axis=0, keepdims=True), columns=cn) perspective_jsd = perspective_jsd_matrix(opinions, nTopics) ps = opinions.keys() for p1, p2 in combinations(ps, 2): op = filter_opinions([p1, p2], opinions) jsd = jsd_for_all_topics(op) idx1 = ps.index(p1) idx2 = ps.index(p2) for t in range(nTopics): yield assert_equal, jsd[t], perspective_jsd[t, idx1, idx2]
def DBSCAN(D, eps, MinPts): noise = [] visited = [] C = [] c_n = -1 for point in D: visited.append(point) # marking point as visited # print point neighbourPts = regionQuery(point, eps, D) if len(neighbourPts) < MinPts: noise.append(point) else: C.append([]) c_n += 1 expandCluster(point, neighbourPts, C, c_n, eps, MinPts, D, visited) print("no. of clusters: ", len(C)) print("length of noise:", len(noise)) for cluster in C: col = [rand(1), rand(1), rand(1)] print(cluster) plt.scatter([i[1] for i in cluster], [i[2] for i in cluster], color=col) plt.show()
def test_fit(): npr.seed(1) N = 10 D = 5 gp = GP(D, burnin=5) inputs = npr.rand(N,D) pending = npr.rand(3,D) W = npr.randn(D,1) vals = inputs.dot(W).flatten() + np.sqrt(1e-3)*npr.randn(N) gp.fit(inputs, vals, pending) assert gp.chain_length == 15 assert all([np.all(p.value != p.initial_value) for p in gp.params.values()]) assert len(gp._cache_list) == 10 assert len(gp._hypers_list) == 10 assert len(gp._fantasy_values_list) == 10
def logprob(self, x): lp = 0.0 for prior in self.priors: lp += prior.logprob(x) return lp # class Binomial(AbstractPrior): # def __init__(self, p, n): # self.p = p # self.n = n # def logprob(self, k): # pos = k # neg = self.n-k # with np.errstate(divide='ignore'): # suppress warnings about log(0) # return np.sum( pos[pos>0]*np.log(self.p[pos>0]) ) + np.sum( neg[neg>0]*np.log(1-self.p[neg>0]) ) # def sample(self, n_samples): # return np.sum(npr.rand(n, n_samples) < p, axis=0) # class Bernoulli(Binomial): # def __init__(self, p): # super(Bernoulli, self).__init__(p, 1)
def add_noise_to_string(self, a_string, amount_of_noise): """Add some artificial spelling mistakes to the string""" if rand() < amount_of_noise * len(a_string): # Replace a character with a random character random_char_position = random_randint(len(a_string)) a_string = a_string[:random_char_position] + random_choice(CHARS[:-1]) + a_string[random_char_position + 1:] if rand() < amount_of_noise * len(a_string): # Delete a character random_char_position = random_randint(len(a_string)) a_string = a_string[:random_char_position] + a_string[random_char_position + 1:] if len(a_string) < MAX_INPUT_LEN and rand() < amount_of_noise * len(a_string): # Add a random character random_char_position = random_randint(len(a_string)) a_string = a_string[:random_char_position] + random_choice(CHARS[:-1]) + a_string[random_char_position:] if rand() < amount_of_noise * len(a_string): # Transpose 2 characters random_char_position = random_randint(len(a_string) - 1) a_string = (a_string[:random_char_position] + a_string[random_char_position + 1] + a_string[random_char_position] + a_string[random_char_position + 2:]) return a_string
def test_boolean(self): a = rand(3, 5, 8) V = rand(5, 8) g1 = randint(0, 5, size=15) g2 = randint(0, 8, size=15) V[g1, g2] = -V[g1, g2] assert_((np.array([a[0][V > 0], a[1][V > 0], a[2][V > 0]]) == a[:, V > 0]).all())
def _generate_data_complex(self, n, m): return randn(n, m) + 1.j * rand(n, m)
def _generate_int_data(self, n, m): return (10 * rand(n, m)).astype(np.int64)
def _generate_int32_data(self, n, m): return (10 * rand(n, m)).astype(np.int32) # Now the real test cases
def test_nd(self): x = 20 * rand(10, 20, 30) out1 = x[:, :, 1:] - x[:, :, :-1] out2 = out1[:, :, 1:] - out1[:, :, :-1] out3 = x[1:, :, :] - x[:-1, :, :] out4 = out3[1:, :, :] - out3[:-1, :, :] assert_array_equal(diff(x), out1) assert_array_equal(diff(x, n=2), out2) assert_array_equal(diff(x, axis=0), out3) assert_array_equal(diff(x, n=2, axis=0), out4)
def test_both(self): a = rand(10) mask = a > 0.5 ac = a.copy() c = extract(mask, a) place(a, mask, 0) place(a, mask, c) assert_array_equal(a, ac)
def test_random(self): x = rand(10) bin = np.linspace(x.min(), x.max(), 10) assert_(np.all(digitize(x, bin) != 0))
def test_right_open_random(self): x = rand(10) bins = np.linspace(x.min(), x.max(), 10) assert_(np.all(digitize(x, bins, True) != 10))
def test_simple(self): # check that unwrap removes jumps greather that 2*pi assert_array_equal(unwrap([1, 1 + 2 * np.pi]), [1, 1]) # check that unwrap maintans continuity assert_(np.all(diff(unwrap(rand(10) * 100)) < np.pi))
def test_normed(self): # Check that the integral of the density equals 1. n = 100 v = rand(n) a, b = histogram(v, normed=True) area = np.sum(a * diff(b)) assert_almost_equal(area, 1) # Check with non-constant bin widths (buggy but backwards # compatible) v = np.arange(10) bins = [0, 1, 5, 9, 10] a, b = histogram(v, bins, normed=True) area = np.sum(a * diff(b)) assert_almost_equal(area, 1)
def test_density(self): # Check that the integral of the density equals 1. n = 100 v = rand(n) a, b = histogram(v, density=True) area = np.sum(a * diff(b)) assert_almost_equal(area, 1) # Check with non-constant bin widths v = np.arange(10) bins = [0, 1, 3, 6, 10] a, b = histogram(v, bins, density=True) assert_array_equal(a, .1) assert_equal(np.sum(a * diff(b)), 1) # Variale bin widths are especially useful to deal with # infinities. v = np.arange(10) bins = [0, 1, 3, 6, np.inf] a, b = histogram(v, bins, density=True) assert_array_equal(a, [.1, .1, .1, 0.]) # Taken from a bug report from N. Becker on the numpy-discussion # mailing list Aug. 6, 2010. counts, dmy = np.histogram( [1, 2, 3, 4], [0.5, 1.5, np.inf], density=True) assert_equal(counts, [.25, 0])
def test_weights(self): v = rand(100) w = np.ones(100) * 5 a, b = histogram(v) na, nb = histogram(v, normed=True) wa, wb = histogram(v, weights=w) nwa, nwb = histogram(v, weights=w, normed=True) assert_array_almost_equal(a * 5, wa) assert_array_almost_equal(na, nwa) # Check weights are properly applied. v = np.linspace(0, 10, 10) w = np.concatenate((np.zeros(5), np.ones(5))) wa, wb = histogram(v, bins=np.arange(11), weights=w) assert_array_almost_equal(wa, w) # Check with integer weights wa, wb = histogram([1, 2, 2, 4], bins=4, weights=[4, 3, 2, 1]) assert_array_equal(wa, [4, 5, 0, 1]) wa, wb = histogram( [1, 2, 2, 4], bins=4, weights=[4, 3, 2, 1], normed=True) assert_array_almost_equal(wa, np.array([4, 5, 0, 1]) / 10. / 3. * 4) # Check weights with non-uniform bin widths a, b = histogram( np.arange(9), [0, 1, 3, 6, 10], weights=[2, 1, 1, 1, 1, 1, 1, 1, 1], density=True) assert_almost_equal(a, [.2, .1, .1, .075])
def test_shape_3d(self): # All possible permutations for bins of different lengths in 3D. bins = ((5, 4, 6), (6, 4, 5), (5, 6, 4), (4, 6, 5), (6, 5, 4), (4, 5, 6)) r = rand(10, 3) for b in bins: H, edges = histogramdd(r, b) assert_(H.shape == b)
def test_shape_4d(self): # All possible permutations for bins of different lengths in 4D. bins = ((7, 4, 5, 6), (4, 5, 7, 6), (5, 6, 4, 7), (7, 6, 5, 4), (5, 7, 6, 4), (4, 6, 7, 5), (6, 5, 7, 4), (7, 5, 4, 6), (7, 4, 6, 5), (6, 4, 7, 5), (6, 7, 5, 4), (4, 6, 5, 7), (4, 7, 5, 6), (5, 4, 6, 7), (5, 7, 4, 6), (6, 7, 4, 5), (6, 5, 4, 7), (4, 7, 6, 5), (4, 5, 6, 7), (7, 6, 4, 5), (5, 4, 7, 6), (5, 6, 7, 4), (6, 4, 5, 7), (7, 5, 6, 4)) r = rand(10, 4) for b in bins: H, edges = histogramdd(r, b) assert_(H.shape == b)
