def PCA(data, num_components=None):
# mean center the data
data -= data.mean(axis=0)
# calculate the covariance matrix
R = np.cov(data, rowvar=False)
# calculate eigenvectors & eigenvalues of the covariance matrix
# use 'eigh' rather than 'eig' since R is symmetric,
# the performance gain is substantial
V, E = np.linalg.eigh(R)
# sort eigenvalue in decreasing order
idx = np.argsort(V)[::-1]
E = E[:,idx]
# sort eigenvectors according to same index
V = V[idx]
# select the first n eigenvectors (n is desired dimension
# of rescaled data array, or dims_rescaled_data)
E = E[:, :num_components]
# carry out the transformation on the data using eigenvectors
# and return the re-scaled data, eigenvalues, and eigenvectors
return np.dot(E.T, data.T).T, V, E
python类dot()的实例源码
def backPropagate(Z1, Z2, y, W2, b2):
## YOUR CODE HERE ##
E2 = 0
E1 = 0
Eb1 = 0
# E2 is the error in output layer. To find it we should exract estimated value from actual output.
# We should find 5 error because there are 5 node in output layer.
E2 = Z2 - y
## E1 is the error in the hidden layer. To find it we should use the error that we found in output layer and the weights between
## output and hidden layer
## We should find 30 error because there are 30 node in hidden layer.
E1 = np.dot(W2, np.transpose(E2))
## Eb1 is the error bias for hidden layer. To find it we should use the error that we found in output layer and the weights between
## output and bias layer
## We should find 1 error because there are 1 bias node in hidden layer.
Eb1 = np.dot(b2, np.transpose(E2))
####################
return E2, E1, Eb1
# calculate the gradients for weights between units and the bias weights
def get_nodal_differentiation_matrix(order,
s2c=None,c2s=None,
Dmodal=None):
"""
Returns the differentiation matrix for the first derivative
in the nodal basis
It goes without saying that this differentiation matrix is for the
reference cell.
"""
if Dmodal is None:
Dmodal = get_modal_differentiation_matrix(order)
if s2c is None or c2s is None:
s2c,c2s = get_vandermonde_matrices(order)
return np.dot(s2c,np.dot(Dmodal,c2s))
# ======================================================================
# Operators Outside Reference Cell
# ======================================================================
def differentiate(self,grid_func,orderx,ordery):
"""Given a grid function defined on the colocation points,
differentiate it up to the appropriate order in each direction.
"""
assert type(orderx) is int
assert type(ordery) is int
assert orderx >= 0
assert ordery >= 0
if orderx > 0:
df = np.dot(self.stencil_x.PD,grid_func)
return self.differentiate(df,orderx-1,ordery)
if ordery > 0:
df = np.dot(grid_func,self.stencil_y.PD.transpose())
return self.differentiate(df,orderx,ordery-1)
#if orderx == 0 and ordery == 0:
return grid_func
def fit(self, graphs, y=None):
rnd = check_random_state(self.random_state)
n_samples = len(graphs)
# get basis vectors
if self.n_components > n_samples:
n_components = n_samples
else:
n_components = self.n_components
n_components = min(n_samples, n_components)
inds = rnd.permutation(n_samples)
basis_inds = inds[:n_components]
basis = []
for ind in basis_inds:
basis.append(graphs[ind])
basis_kernel = self.kernel(basis, basis, **self._get_kernel_params())
# sqrt of kernel matrix on basis vectors
U, S, V = svd(basis_kernel)
S = np.maximum(S, 1e-12)
self.normalization_ = np.dot(U * 1. / np.sqrt(S), V)
self.components_ = basis
self.component_indices_ = inds
return self
def rhoA(self):
# rhoA
rhoA = pd.DataFrame(0, index=np.arange(1), columns=self.latent)
for i in range(self.lenlatent):
weights = pd.DataFrame(self.outer_weights[self.latent[i]])
weights = weights[(weights.T != 0).any()]
result = pd.DataFrame.dot(weights.T, weights)
result_ = pd.DataFrame.dot(weights, weights.T)
S = self.data_[self.Variables['measurement'][
self.Variables['latent'] == self.latent[i]]]
S = pd.DataFrame.dot(S.T, S) / S.shape[0]
numerador = (
np.dot(np.dot(weights.T, (S - np.diag(np.diag(S)))), weights))
denominador = (
(np.dot(np.dot(weights.T, (result_ - np.diag(np.diag(result_)))), weights)))
rhoA_ = ((result)**2) * (numerador / denominador)
if(np.isnan(rhoA_.values)):
rhoA[self.latent[i]] = 1
else:
rhoA[self.latent[i]] = rhoA_.values
return rhoA.T
def _ikf_iteration(self, x, n, ranges, h, H, z, estimate, R):
"""Update tracker based on a multi-range message.
Args:
multi_range_msg (uwb.msg.UWBMultiRangeWithOffsets): ROS multi-range message.
Returns:
new_estimate (StateEstimate): Updated position estimate.
"""
new_position = n[0:3]
self._compute_measurements_and_jacobians(ranges, new_position, h, H, z)
res = z - h
S = np.dot(np.dot(H, estimate.covariance), H.T) + R
K = np.dot(estimate.covariance, self._solve_equation_least_squares(S.T, H).T)
mahalanobis = np.sqrt(np.dot(self._solve_equation_least_squares(S.T, res).T, res))
if res.size not in self.outlier_thresholds:
self.outlier_thresholds[res.size] = scipy.stats.chi2.isf(self.outlier_threshold_quantile, res.size)
outlier_threshold = self.outlier_thresholds[res.size]
if mahalanobis < outlier_threshold:
n = x + np.dot(K, (res - np.dot(H, x - n)))
outlier_flag = False
else:
outlier_flag = True
return n, K, outlier_flag
def normalized_distance(_a, _b):
"""Compute normalized distance between two points.
Computes 1 - a * b / ( ||a|| * ||b||)
Args:
_a (numpy.ndarray): array of size m
_b (numpy.ndarray): array of size m
Returns:
normalized distance between signatures (float)
Examples:
>>> a = gis.generate_signature('https://upload.wikimedia.org/wikipedia/commons/thumb/e/ec/Mona_Lisa,_by_Leonardo_da_Vinci,_from_C2RMF_retouched.jpg/687px-Mona_Lisa,_by_Leonardo_da_Vinci,_from_C2RMF_retouched.jpg')
>>> b = gis.generate_signature('https://pixabay.com/static/uploads/photo/2012/11/28/08/56/mona-lisa-67506_960_720.jpg')
>>> gis.normalized_distance(a, b)
0.0332806110382
"""
# return (1.0 - np.dot(_a, _b) / (np.linalg.norm(_a) * np.linalg.norm(_b)))
return np.dot(_a, _b) / (np.linalg.norm(_a) * np.linalg.norm(_b))
def observed_perplexity(self, counts):
"""Compute perplexity = exp(entropy) of observed variables.
Perplexity is an information theoretic measure of the number of
clusters or latent classes. Perplexity is a real number in the range
[1, M], where M is model_num_clusters.
Args:
counts: A [V]-shaped array of multinomial counts.
Returns:
A [V]-shaped numpy array of perplexity.
"""
V, E, M, R = self._VEMR
if counts is not None:
counts = np.ones(V, dtype=np.int8)
assert counts.shape == (V, )
assert counts.dtype == np.int8
assert np.all(counts > 0)
observed_entropy = np.empty(V, dtype=np.float32)
for v in range(V):
beg, end = self._ragged_index[v:v + 2]
probs = np.dot(self._feat_cond[beg:end, :], self._vert_probs[v, :])
observed_entropy[v] = multinomial_entropy(probs, counts[v])
return np.exp(observed_entropy)
def get_covariance(self):
"""Compute data covariance with the generative model.
``cov = components_.T * S**2 * components_ + sigma2 * eye(n_features)``
where S**2 contains the explained variances.
Returns
-------
cov : array, shape=(n_features, n_features)
Estimated covariance of data.
"""
components_ = self.components_
exp_var = self.explained_variance_
if self.whiten:
components_ = components_ * np.sqrt(exp_var[:, np.newaxis])
exp_var_diff = np.maximum(exp_var - self.noise_variance_, 0.)
cov = np.dot(components_.T * exp_var_diff, components_)
cov.flat[::len(cov) + 1] += self.noise_variance_ # modify diag inplace
return cov
def transform(self, X):
"""Apply the dimensionality reduction on X.
X is projected on the first principal components previous extracted
from a training set.
Parameters
----------
X : array-like, shape (n_samples, n_features)
New data, where n_samples is the number of samples
and n_features is the number of features.
Returns
-------
X_new : array-like, shape (n_samples, n_components)
"""
check_is_fitted(self, 'mean_')
X = check_array(X)
if self.mean_ is not None:
X = X - self.mean_
X_transformed = np.dot(X, self.components_.T)
if self.whiten:
X_transformed /= np.sqrt(self.explained_variance_)
return X_transformed
def score_samples(self, X):
"""Return the log-likelihood of each sample
See. "Pattern Recognition and Machine Learning"
by C. Bishop, 12.2.1 p. 574
or http://www.miketipping.com/papers/met-mppca.pdf
Parameters
----------
X: array, shape(n_samples, n_features)
The data.
Returns
-------
ll: array, shape (n_samples,)
Log-likelihood of each sample under the current model
"""
check_is_fitted(self, 'mean_')
X = check_array(X)
Xr = X - self.mean_
n_features = X.shape[1]
log_like = np.zeros(X.shape[0])
precision = self.get_precision()
log_like = -.5 * (Xr * (np.dot(Xr, precision))).sum(axis=1)
log_like -= .5 * (n_features * log(2. * np.pi)
- fast_logdet(precision))
return log_like
def optSigma2(self, U, s, y, covars, logdetXX, reml, ldeltamin=-5, ldeltamax=5):
#Prepare required matrices
Uy = U.T.dot(y).flatten()
UX = U.T.dot(covars)
if (U.shape[1] < U.shape[0]):
UUX = covars - U.dot(UX)
UUy = y - U.dot(Uy)
UUXUUX = UUX.T.dot(UUX)
UUXUUy = UUX.T.dot(UUy)
UUyUUy = UUy.T.dot(UUy)
else: UUXUUX, UUXUUy, UUyUUy = None, None, None
numIndividuals = U.shape[0]
ldeltaopt_glob = optimize.minimize_scalar(self.negLLevalLong, bounds=(-5, 5), method='Bounded', args=(s, Uy, UX, logdetXX, UUXUUX, UUXUUy, UUyUUy, numIndividuals, reml)).x
ll, sig2g, beta, r2 = self.negLLevalLong(ldeltaopt_glob, s, Uy, UX, logdetXX, UUXUUX, UUXUUy, UUyUUy, numIndividuals, reml, returnAllParams=True)
sig2e = np.exp(ldeltaopt_glob) * sig2g
return sig2g, sig2e, beta, ll
def BorthogonalityTest(B, U):
"""
Test the frobenious norm of U^TBU - I_k
"""
BU = np.zeros(U.shape)
Bu = Vector()
u = Vector()
B.init_vector(Bu,0)
B.init_vector(u,1)
nvec = U.shape[1]
for i in range(0,nvec):
u.set_local(U[:,i])
B.mult(u,Bu)
BU[:,i] = Bu.get_local()
UtBU = np.dot(U.T, BU)
err = UtBU - np.eye(nvec, dtype=UtBU.dtype)
print("|| UtBU - I ||_F = ", np.linalg.norm(err, 'fro') )
def _ireduce_linalg(arrays, func, **kwargs):
"""
Yield the cumulative reduction of a linag algebra function
"""
arrays = iter(arrays)
first = next(arrays)
second = next(arrays)
func = partial(func, **kwargs)
accumulator = func(first, second)
yield accumulator
for array in arrays:
# For some reason, np.dot(..., out = accumulator) did not produce results
# that were equal to numpy.linalg.multi_dot
func(accumulator, array, out = accumulator)
yield accumulator
def idot(arrays):
"""
Yields the cumulative array inner product (dot product) of arrays.
Parameters
----------
arrays : iterable
Arrays to be reduced.
Yields
------
online_dot : ndarray
See Also
--------
numpy.linalg.multi_dot : Compute the dot product of two or more arrays in a single function call,
while automatically selecting the fastest evaluation order.
"""
yield from _ireduce_linalg(arrays, np.dot)
def itensordot(arrays, axes = 2):
"""
Yields the cumulative array inner product (dot product) of arrays.
Parameters
----------
arrays : iterable
Arrays to be reduced.
axes : int or (2,) array_like
* integer_like: If an int N, sum over the last N axes of a
and the first N axes of b in order. The sizes of the corresponding axes must match.
* (2,) array_like: Or, a list of axes to be summed over, first sequence applying to a,
second to b. Both elements array_like must be of the same length.
Yields
------
online_tensordot : ndarray
See Also
--------
numpy.tensordot : Compute the tensordot on two tensors.
"""
yield from _ireduce_linalg(arrays, np.tensordot, axes = axes)
def _convert(matrix, arr):
"""Do the color space conversion.
Parameters
----------
matrix : array_like
The 3x3 matrix to use.
arr : array_like
The input array.
Returns
-------
out : ndarray, dtype=float
The converted array.
"""
arr = _prepare_colorarray(arr)
arr = np.swapaxes(arr, 0, -1)
oldshape = arr.shape
arr = np.reshape(arr, (3, -1))
out = np.dot(matrix, arr)
out.shape = oldshape
out = np.swapaxes(out, -1, 0)
return np.ascontiguousarray(out)
def rotate_point_cloud(batch_data):
""" Randomly rotate the point clouds to augument the dataset
rotation is per shape based along up direction
Input:
BxNx3 array, original batch of point clouds
Return:
BxNx3 array, rotated batch of point clouds
"""
rotated_data = np.zeros(batch_data.shape, dtype=np.float32)
for k in range(batch_data.shape[0]):
rotation_angle = np.random.uniform() * 2 * np.pi
cosval = np.cos(rotation_angle)
sinval = np.sin(rotation_angle)
rotation_matrix = np.array([[cosval, 0, sinval],
[0, 1, 0],
[-sinval, 0, cosval]])
shape_pc = batch_data[k, ...]
rotated_data[k, ...] = np.dot(shape_pc.reshape((-1, 3)), rotation_matrix)
return rotated_data
def rotate_point_cloud_by_angle(batch_data, rotation_angle):
""" Rotate the point cloud along up direction with certain angle.
Input:
BxNx3 array, original batch of point clouds
Return:
BxNx3 array, rotated batch of point clouds
"""
rotated_data = np.zeros(batch_data.shape, dtype=np.float32)
for k in range(batch_data.shape[0]):
#rotation_angle = np.random.uniform() * 2 * np.pi
cosval = np.cos(rotation_angle)
sinval = np.sin(rotation_angle)
rotation_matrix = np.array([[cosval, 0, sinval],
[0, 1, 0],
[-sinval, 0, cosval]])
shape_pc = batch_data[k, ...]
rotated_data[k, ...] = np.dot(shape_pc.reshape((-1, 3)), rotation_matrix)
return rotated_data
def computeStep(X, y, theta):
'''YOUR CODE HERE'''
function_result = np.array([0,0], dtype= np.float)
m = float(len(X))
d1 = 0.0
d2 = 0.0
for i in range(len(X)):
h1 = np.dot(theta.transpose(), X[i])
c1 = h1 - y[i]
d1 = d1 + c1
j1 = d1/m
for u in range(len(X)):
h2 = np.dot(theta.transpose(), X[u])
c2 = (h2 - y[u]) * X[u][1]
d2 = d2 + c2
j2 = d2/m
function_result[0] = j1
function_result[1] = j2
return function_result
# Part 4: Implement the cost function calculation
def computeCost(X, y, theta):
'''YOUR CODE HERE'''
m = float(len(X))
d = 0
for i in range(len(X)):
h = np.dot(theta.transpose(), X[i])
c = (h - y[i])
c = (c **2)
d = (d + c)
j = (1.0 / (2 * m)) * d
return j
# Part 5: Prepare the data so that the input X has two columns: first a column of ones to accomodate theta0 and then a column of city population data
def forwardPropagate(X, W1, b1, W2, b2):
## YOUR CODE HERE ##
Z1 = 0
Z2 = 0
S1 = 0
S2 = 0
## Here we should find 2 result: first for input and hidden layer then for hidden and output layer.
## First I found the result for every node in hidden layer then put them into activation function.
S1 = np.dot(X, W1)+ b1
Z1 = calcActivation(S1)
## Second I found the result for every node in output layer then put them into activation function.
S2 = np.dot(Z1, W2) + b2
Z2 = calcActivation(S2)
####################
return Z1, Z2
# calculate the cost
def calcGrads(X, Z1, Z2, E1, E2, Eb1):
## YOUR CODE HERE ##
d_W1 = 0
d_b1 = 0
d_W2 = 0
d_b2 = 0
## In here we should the derivatives for gradients. To find derivative, we should multiply.
# d_w2 is the derivative for weights between hidden layer and the output layer.
d_W2 = np.dot(np.transpose(E2), Z1)
# d_w1 is the derivative for weights between hidden layer and the input layer.
d_W1 = np.dot(E1, X)
# d_b2 is the derivative for weights between hidden layer bias and the output layer.
d_b2 = np.dot(np.transpose(E2), Eb1)
# d_b1 is the derivative for weights between hidden layer bias and the input layer.
d_b1 = np.dot(np.transpose(E1), 1)
####################
return d_W1, d_W2, d_b1, d_b2
# update the weights between units and the bias weights using a learning rate of alpha
def doesnt_match(self, words):
"""
Which word from the given list doesn't go with the others?
Example::
>>> trained_model.doesnt_match("breakfast cereal dinner lunch".split())
'cereal'
"""
words = [word for word in words if word in self.vocab] # filter out OOV words
logger.debug("using words %s" % words)
if not words:
raise ValueError("cannot select a word from an empty list")
# which word vector representation is furthest away from the mean?
selection = self.syn0norm[[self.vocab[word].index for word in words]]
mean = np.mean(selection, axis=0)
sim = np.dot(selection, mean / np.linalg.norm(mean))
return words[np.argmin(sim)]
def __mul__(self, other):
"""
Left-multiply RigidTransform with another rigid transform
Two variants:
RigidTransform: Identical to oplus operation
ndarray: transform [N x 3] point set (X_2 = p_21 * X_1)
"""
if isinstance(other, DualQuaternion):
return DualQuaternion.from_dq(other.real * self.real,
other.dual * self.real + other.real * self.dual)
elif isinstance(other, float):
return DualQuaternion.from_dq(self.real * other, self.dual * other)
# elif isinstance(other, nd.array):
# X = np.hstack([other, np.ones((len(other),1))]).T
# return (np.dot(self.matrix, X).T)[:,:3]
else:
raise TypeError('__mul__ typeerror {:}'.format(type(other)))
def quaternion_matrix(quaternion):
"""Return homogeneous rotation matrix from quaternion.
>>> R = quaternion_matrix([0.06146124, 0, 0, 0.99810947])
>>> numpy.allclose(R, rotation_matrix(0.123, (1, 0, 0)))
True
"""
q = numpy.array(quaternion[:4], dtype=numpy.float64, copy=True)
nq = numpy.dot(q, q)
if nq < _EPS:
return numpy.identity(4)
q *= math.sqrt(2.0 / nq)
q = numpy.outer(q, q)
return numpy.array((
(1.0-q[1, 1]-q[2, 2], q[0, 1]-q[2, 3], q[0, 2]+q[1, 3], 0.0),
( q[0, 1]+q[2, 3], 1.0-q[0, 0]-q[2, 2], q[1, 2]-q[0, 3], 0.0),
( q[0, 2]-q[1, 3], q[1, 2]+q[0, 3], 1.0-q[0, 0]-q[1, 1], 0.0),
( 0.0, 0.0, 0.0, 1.0)
), dtype=numpy.float64)
def sampson_error(F, pts1, pts2):
"""
Computes the sampson error for F, and points pts1, pts2. Sampson
error is the first order approximation to the geometric error.
Remember that this is a squared error.
(x'^{T} * F * x)^2
-----------------
(F * x)_1^2 + (F * x)_2^2 + (F^T * x')_1^2 + (F^T * x')_2^2
where (F * x)_i^2 is the square of the i-th entry of the vector Fx
"""
x1, x2 = unproject_points(pts1).T, unproject_points(pts2).T
Fx1 = np.dot(F, x1)
Fx2 = np.dot(F, x2)
# Sampson distance as error measure
denom = Fx1[0]**2 + Fx1[1]**2 + Fx2[0]**2 + Fx2[1]**2
return ( np.diag(x1.T.dot(Fx2)) )**2 / denom
def getMedianDistanceBetweenSamples(self, sampleSet=None) :
"""
Jaakkola's heuristic method for setting the width parameter of the Gaussian
radial basis function kernel is to pick a quantile (usually the median) of
the distribution of Euclidean distances between points having different
labels.
Reference:
Jaakkola, M. Diekhaus, and D. Haussler. Using the Fisher kernel method to detect
remote protein homologies. In T. Lengauer, R. Schneider, P. Bork, D. Brutlad, J.
Glasgow, H.- W. Mewes, and R. Zimmer, editors, Proceedings of the Seventh
International Conference on Intelligent Systems for Molecular Biology.
"""
numrows = sampleSet.shape[0]
samples = sampleSet
G = sum((samples * samples), 1)
Q = numpy.tile(G[:, None], (1, numrows))
R = numpy.tile(G, (numrows, 1))
distances = Q + R - 2 * numpy.dot(samples, samples.T)
distances = distances - numpy.tril(distances)
distances = distances.reshape(numrows**2, 1, order="F").copy()
return numpy.sqrt(0.5 * numpy.median(distances[distances > 0]))
def refit_model(self):
"""Learns a new surrogate model using the data observed so far.
"""
# only fit the model if there is data for it.
if len(self.known_models) > 0:
self._build_feature_maps(self.known_models, self.ngram_maxlen, self.thres)
X = sp.vstack([ self._compute_features(mdl)
for mdl in self.known_models], "csr")
y = np.array(self.known_scores, dtype='float64')
#A = np.dot(X.T, X) + lamb * np.eye(X.shape[1])
#b = np.dot(X.T, y)
self.surr_model = lm.Ridge(self.lamb_ridge)
self.surr_model.fit(X, y)
# NOTE: if the search space has holes, it break. needs try/except module.
