def array2PIL(arr, size):
mode = 'RGBA'
arr = arr.reshape(arr.shape[0]*arr.shape[1], arr.shape[2])
if len(arr[0]) == 3:
arr = numpy.c_[arr, 255*numpy.ones((len(arr),1), numpy.uint8)]
return Image.frombuffer(mode, size, arr.tostring(), 'raw', mode, 0, 1)
python类c_()的实例源码
def best_fit_plane(self):
"""Fits a plane to the point cloud using least squares.
Returns
-------
:obj:`tuple` of :obj:`numpy.ndarray` of float
A normal vector to and point in the fitted plane.
"""
X = np.c_[self.x_coords, self.y_coords, np.ones(self.num_points)]
y = self.z_coords
A = X.T.dot(X)
b = X.T.dot(y)
w = np.linalg.inv(A).dot(b)
n = np.array([w[0], w[1], -1])
n = n / np.linalg.norm(n)
n = Direction(n, self._frame)
x0 = self.mean()
return n, x0
def predict_bbox_regressor(model, feat, ex_boxes):
if ex_boxes.size == 0:
return np.array([]).reshape(-1, 4)
# predict regression targets
Y = np.dot(feat, model.Beta[:-1]) + model.Beta[-1]
# invert transformation
Y = dot(Y, model.T_inv)
# read out prediction
dst_size = Y[:, 2:]
dst_ctr = Y[:, 2:]
src_size = ex_boxes[:, 2:]
src_ctr = ex_boxes[:, :2] + 0.5 * src_size
pred_size = np.exp(dst_size) * src_size
pred_ctr = dst_ctr * src_ctr + src_ctr
pred = np.c_[pred_ctr - 0.5 * pred_size, pred_size]
return pred
def get_examples(bbox, gt):
# compute overlap ratio
O = overlap_ratio(bbox, gt)
# compute answer
src_size = bbox[:, 2:]
src_ctr = bbox[:, :2] + 0.5 * src_size
gt_size = gt[2:]
gt_ctr = gt[:2] + 0.5 * gt_size
dst_size = np.log(gt_size / src_size)
dst_ctr = (gt_ctr - src_ctr) * 1.0 / src_ctr
Y = np.c_[dst_ctr, dst_size]
return Y, O
def computeRotMatrix(self,Phi=False):
#######################################
# COMPUTE ROTATION MATRIX SUCH THAT m(t) = A*L(t)*A'*Hp
# Default set such that phi1,phi2 = 0 is UXO pointed towards North
if Phi is False:
phi1 = np.radians(self.phi[0])
phi2 = np.radians(self.phi[1])
phi3 = np.radians(self.phi[2])
else:
phi1 = np.radians(Phi[0]) # Roll (CCW)
phi2 = np.radians(Phi[1]) # Inclination (+ve is nose pointing down)
phi3 = np.radians(Phi[2]) # Declination (degrees CW from North)
# A1 = np.r_[np.c_[np.cos(phi1),-np.sin(phi1),0.],np.c_[np.sin(phi1),np.cos(phi1),0.],np.c_[0.,0.,1.]] # CCW Rotation about z-axis
# A2 = np.r_[np.c_[1.,0.,0.],np.c_[0.,np.cos(phi2),np.sin(phi2)],np.c_[0.,-np.sin(phi2),np.cos(phi2)]] # CW Rotation about x-axis (rotates towards North)
# A3 = np.r_[np.c_[np.cos(phi3),-np.sin(phi3),0.],np.c_[np.sin(phi3),np.cos(phi3),0.],np.c_[0.,0.,1.]] # CCW Rotation about z-axis (direction of head of object)
A1 = np.r_[np.c_[np.cos(phi1),np.sin(phi1),0.],np.c_[-np.sin(phi1),np.cos(phi1),0.],np.c_[0.,0.,1.]] # CW Rotation about z-axis
A2 = np.r_[np.c_[1.,0.,0.],np.c_[0.,np.cos(phi2),np.sin(phi2)],np.c_[0.,-np.sin(phi2),np.cos(phi2)]] # CW Rotation about x-axis (rotates towards North)
A3 = np.r_[np.c_[np.cos(phi3),np.sin(phi3),0.],np.c_[-np.sin(phi3),np.cos(phi3),0.],np.c_[0.,0.,1.]] # CW Rotation about z-axis (direction of head of object)
return np.dot(A3,np.dot(A2,A1))
def defineSensorLoc(self,XYZ):
#############################################
# DEFINE TRANSMITTER AND RECEIVER LOCATIONS
#
# XYZ: N X 3 array containing transmitter center locations
# **NOTE** for this sensor, we know where the receivers are relative to transmitters
self.TxLoc = XYZ
dx,dy = np.meshgrid([-0.8,-0.4,0.,0.4,0.8],[-0.8,-0.4,0.,0.4,0.8])
dx = mkvc(dx)
dy = mkvc(dy)
N = np.shape(XYZ)[0]
X = np.kron(XYZ[:,0],np.ones((25))) + np.kron(np.ones((N)),dx)
Y = np.kron(XYZ[:,1],np.ones((25))) + np.kron(np.ones((N)),dy)
Z = np.kron(XYZ[:,2],np.ones((25)))
self.RxLoc = np.c_[X,Y,Z]
self.TxID = np.kron(np.arange(1,np.shape(XYZ)[0]+1),np.ones((25)))
self.RxComp = np.kron(3*np.ones(np.shape(XYZ)[0]),np.ones((25)))
def updatePolarizations(self,r0,UB):
# Set operator and solution array
Hp = self.computeHp(r0=r0)
Brx = self.computeBrx(r0=r0)
P = self.computeP(Hp,Brx)
dunc = self.dunc
dobs = self.dobs
K = np.shape(dobs)[1]
q = np.zeros((6,K))
lb = np.zeros(6)
ub = UB*np.ones(6)
for kk in range(0,K):
LHS = P/np.c_[dunc[:,kk],dunc[:,kk],dunc[:,kk],dunc[:,kk],dunc[:,kk],dunc[:,kk]]
RHS = dobs[:,kk]/dunc[:,kk]
Sol = op.lsq_linear(LHS,RHS,bounds=(lb,ub),tol=1e-5)
q[:,kk] = Sol.x
self.q = q
def updatePolarizations(self,r0,UB):
# Set operator and solution array
Hp = self.computeHp(r0=r0)
Brx = self.computeBrx(r0=r0)
P = self.computeP(Hp,Brx)
dunc = self.dunc
dobs = self.dobs
K = np.shape(dobs)[1]
q = np.zeros((6,K))
lb = np.zeros(6)
ub = UB*np.ones(6)
for kk in range(0,K):
LHS = P/np.c_[dunc[:,kk],dunc[:,kk],dunc[:,kk],dunc[:,kk],dunc[:,kk],dunc[:,kk]]
RHS = dobs[:,kk]/dunc[:,kk]
Sol = op.lsq_linear(LHS,RHS,bounds=(lb,ub),tol=1e-7)
q[:,kk] = Sol.x
self.q = q
def analytic_infinite_wire(obsloc,wireloc,orientation,I=1.):
"""
Compute the response of an infinite wire with orientation 'orientation'
and current I at the obsvervation locations obsloc
Output:
B: magnetic field [Bx,By,Bz]
"""
n,d = obsloc.shape
t,d = wireloc.shape
d = np.sqrt(np.dot(obsloc**2.,np.ones([d,t]))+np.dot(np.ones([n,d]),(wireloc.T)**2.)
- 2.*np.dot(obsloc,wireloc.T))
distr = np.amin(d, axis=1, keepdims = True)
idxmind = d.argmin(axis=1)
r = obsloc - wireloc[idxmind]
orient = np.c_[[orientation for i in range(obsloc.shape[0])]]
B = (mu_0*I)/(2*np.pi*(distr**2.))*np.cross(orientation,r)
return B
def predict_with_glm(X, y, model):
""" Predict number of mutation with GLM.
Args:
X (np.array): feature matrix.
y (pd.df): response.
model (dict): model meta-data.
Returns:
np.array: array of predictions.
"""
# Add const. to X
X = np.c_[X, np.ones(X.shape[0])]
if model['model_name'] == 'Binomial':
pred = np.array(model['model'].predict(X) * y.length * y.N)
elif model['model_name'] == 'NegativeBinomial':
pred = np.array(model['model'].predict(X, exposure=(y.length * y.N).values + 1))
else:
sys.stderr.write('Wrong model name in model info: {}. Need Binomial or NegativeBinomial.'.format(model['model_name']))
sys.exit(1)
return pred
def normals(self):
"""Face Normals
:rtype: numpy.array
:return: normals, (sum(nF), dim)
"""
if self.dim == 2:
nX = np.c_[
np.ones(self.nFx), np.zeros(self.nFx)
]
nY = np.c_[
np.zeros(self.nFy), np.ones(self.nFy)
]
return np.r_[nX, nY]
elif self.dim == 3:
nX = np.c_[
np.ones(self.nFx), np.zeros(self.nFx), np.zeros(self.nFx)
]
nY = np.c_[
np.zeros(self.nFy), np.ones(self.nFy), np.zeros(self.nFy)
]
nZ = np.c_[
np.zeros(self.nFz), np.zeros(self.nFz), np.ones(self.nFz)
]
return np.r_[nX, nY, nZ]
def tangents(self):
"""Edge Tangents
:rtype: numpy.array
:return: normals, (sum(nE), dim)
"""
if self.dim == 2:
tX = np.c_[
np.ones(self.nEx), np.zeros(self.nEx)
]
tY = np.c_[
np.zeros(self.nEy), np.ones(self.nEy)
]
return np.r_[tX, tY]
elif self.dim == 3:
tX = np.c_[
np.ones(self.nEx), np.zeros(self.nEx), np.zeros(self.nEx)
]
tY = np.c_[
np.zeros(self.nEy), np.ones(self.nEy), np.zeros(self.nEy)
]
tZ = np.c_[
np.zeros(self.nEz), np.zeros(self.nEz), np.ones(self.nEz)
]
return np.r_[tX, tY, tZ]
def test_invPropertyTensor2D(self):
M = discretize.TensorMesh([6, 6])
a1 = np.random.rand(M.nC)
a2 = np.random.rand(M.nC)
a3 = np.random.rand(M.nC)
prop1 = a1
prop2 = np.c_[a1, a2]
prop3 = np.c_[a1, a2, a3]
for prop in [4, prop1, prop2, prop3]:
b = invPropertyTensor(M, prop)
A = makePropertyTensor(M, prop)
B1 = makePropertyTensor(M, b)
B2 = invPropertyTensor(M, prop, returnMatrix=True)
Z = B1*A - sp.identity(M.nC*2)
self.assertTrue(np.linalg.norm(Z.todense().ravel(), 2) < TOL)
Z = B2*A - sp.identity(M.nC*2)
self.assertTrue(np.linalg.norm(Z.todense().ravel(), 2) < TOL)
def test_TensorType3D(self):
M = discretize.TensorMesh([6, 6, 7])
a1 = np.random.rand(M.nC)
a2 = np.random.rand(M.nC)
a3 = np.random.rand(M.nC)
a4 = np.random.rand(M.nC)
a5 = np.random.rand(M.nC)
a6 = np.random.rand(M.nC)
prop1 = a1
prop2 = np.c_[a1, a2, a3]
prop3 = np.c_[a1, a2, a3, a4, a5, a6]
for ii, prop in enumerate([4, prop1, prop2, prop3]):
self.assertTrue(TensorType(M, prop) == ii)
self.assertRaises(Exception, TensorType, M, np.c_[a1, a2, a3, a3])
self.assertTrue(TensorType(M, None) == -1)
def test_invPropertyTensor3D(self):
M = discretize.TensorMesh([6, 6, 6])
a1 = np.random.rand(M.nC)
a2 = np.random.rand(M.nC)
a3 = np.random.rand(M.nC)
a4 = np.random.rand(M.nC)
a5 = np.random.rand(M.nC)
a6 = np.random.rand(M.nC)
prop1 = a1
prop2 = np.c_[a1, a2, a3]
prop3 = np.c_[a1, a2, a3, a4, a5, a6]
for prop in [4, prop1, prop2, prop3]:
b = invPropertyTensor(M, prop)
A = makePropertyTensor(M, prop)
B1 = makePropertyTensor(M, b)
B2 = invPropertyTensor(M, prop, returnMatrix=True)
Z = B1*A - sp.identity(M.nC*3)
self.assertTrue(np.linalg.norm(Z.todense().ravel(), 2) < TOL)
Z = B2*A - sp.identity(M.nC*3)
self.assertTrue(np.linalg.norm(Z.todense().ravel(), 2) < TOL)
def getError(self):
funR = lambda r, z: np.sin(2.*np.pi*r)
funZ = lambda r, z: np.sin(2.*np.pi*z)
sol = lambda r, t, z: (2*np.pi*r*np.cos(2*np.pi*r) + np.sin(2*np.pi*r))/r + 2*np.pi*np.cos(2*np.pi*z)
Fc = cylF2(self.M, funR, funZ)
Fc = np.c_[Fc[:, 0], np.zeros(self.M.nF), Fc[:, 1]]
F = self.M.projectFaceVector(Fc)
divF = self.M.faceDiv.dot(F)
divF_ana = call3(sol, self.M.gridCC)
err = np.linalg.norm((divF-divF_ana), np.inf)
return err
def getError(self):
funR = lambda r, z: np.sin(2.*np.pi*z) * np.cos(np.pi*r)
funZ = lambda r, z: np.sin(3.*np.pi*z) * np.cos(2.*np.pi*r)
Fc = cylF2(self.M, funR, funZ)
Fc = np.c_[Fc[:, 0], np.zeros(self.M.nF), Fc[:, 1]]
F = self.M.projectFaceVector(Fc)
aveF = self.M.aveF2CCV * F
aveF_anaR = funR(self.M.gridCC[:, 0], self.M.gridCC[:, 2])
aveF_anaZ = funZ(self.M.gridCC[:, 0], self.M.gridCC[:, 2])
aveF_ana = np.hstack([aveF_anaR, aveF_anaZ])
err = np.linalg.norm((aveF-aveF_ana), np.inf)
return err
def surface_points(self, grid_basis=True):
"""Returns the points on the surface.
Parameters
----------
grid_basis : bool
If False, the surface points are transformed to the world frame.
If True (default), the surface points are left in grid coordinates.
Returns
-------
:obj:`tuple` of :obj:`numpy.ndarray` of int, :obj:`numpy.ndarray` of float
The points on the surface and the signed distances at those points.
"""
surface_points = np.where(np.abs(self.data_) < self.surface_thresh_)
x = surface_points[0]
y = surface_points[1]
z = surface_points[2]
surface_points = np.c_[x, np.c_[y, z]]
surface_vals = self.data_[surface_points[:,0], surface_points[:,1], surface_points[:,2]]
if not grid_basis:
surface_points = self.transform_pt_grid_to_obj(surface_points.T)
surface_points = surface_points.T
return surface_points, surface_vals
def _check_freq(f):
"""Check the frequency definition."""
f = np.atleast_2d(np.asarray(f))
#
if len(f.reshape(-1)) == 1:
raise ValueError("The length of f should at least be 2.")
elif 2 in f.shape: # f of shape (N, 2) or (2, N)
if f.shape[1] is not 2:
f = f.T
elif np.squeeze(f).shape == (4,): # (fstart, fend, fwidth, fstep)
f = _pair_vectors(*tuple(np.squeeze(f)))
else: # Sequential
f = f.reshape(-1)
f.sort()
f = np.c_[f[0:-1], f[1::]]
return f
def convert_mask_to_locations(mask):
""" Return the converted Cartesian mask as sampling locations.
Parameters
----------
mask: np.ndarray, {0,1}
2D matrix, not necessarly a square matrix.
Returns
-------
samples_locations: np.ndarray
list of the samples between [-0.5, 0.5[.
"""
row, col = np.where(mask == 1)
row = row.astype("float") / mask.shape[0] - 0.5
col = col.astype("float") / mask.shape[1] - 0.5
return np.c_[row, col]
def _test_double_optimization():
"""Test double optimization on a simple example."""
# A simple sparse-sum function
X = [[1, 2], [3, 4], [5, 6]]
y = [sum(x) for x in X]
T = [[7, 8], [9, 10], [2, 1]]
# noisy variables
np.random.seed(0)
X = np.c_[X, np.random.random((3, 100))]
T = np.c_[T, np.random.random((3, 100))]
# Select the first 2 variables and calculate a linear model on them
dstep = DoubleStepEstimator(Lasso(tau=1.0), RidgeRegression(mu=0.0)).train(X, y)
# Coefficients
lasso = dstep.selector
ridge = dstep.estimator
assert_array_almost_equal([0.90635646, 0.90635646], lasso.beta[:2])
assert_array_almost_equal([1.0, 1.0], ridge.beta)
assert_array_almost_equal([1.0, 1.0], dstep.beta[:2])
# Prediction
y_ = dstep.predict(T)
assert_array_almost_equal([15., 19., 3.], y_)
def generate_hills(width, height, nhills):
'''
@param width float, terrain width
@param height float, terrain height
@param nhills int, #hills to gen. #hills actually generted is sqrt(nhills)^2
'''
# setup coordinate grid
xmin, xmax = -width/2.0, width/2.0
ymin, ymax = -height/2.0, height/2.0
x, y = np.mgrid[xmin:xmax:STEP, ymin:ymax:STEP]
pos = np.empty(x.shape + (2,))
pos[:, :, 0] = x; pos[:, :, 1] = y
# generate hilltops
xm, ym = np.mgrid[xmin:xmax:width/np.sqrt(nhills), ymin:ymax:height/np.sqrt(nhills)]
mu = np.c_[xm.flat, ym.flat]
sigma = float(width*height)/(nhills*8)
for i in range(mu.shape[0]):
mu[i] = multivariate_normal.rvs(mean=mu[i], cov=sigma)
# generate hills
sigma = sigma + sigma*np.random.rand(mu.shape[0])
rvs = [ multivariate_normal(mu[i,:], cov=sigma[i]) for i in range(mu.shape[0]) ]
hfield = np.max([ rv.pdf(pos) for rv in rvs ], axis=0)
return x, y, hfield
def generate_data(sample_size=200, pd=[[0.4, 0.4], [0.1, 0.1]]):
pd = np.array(pd)
pd /= pd.sum()
offset = 50
bins = np.r_[np.zeros((1,)), np.cumsum(pd)]
bin_counts = np.histogram(np.random.rand(sample_size), bins)[0]
data = np.empty((0, 2))
targets = []
for ((i, j), p), count in zip(np.ndenumerate(pd), bin_counts):
xs = np.random.uniform(low=0.0, high=50.0, size=count) + j * offset
ys = np.random.uniform(low=0.0, high=50.0, size=count) + -i * offset
data = np.vstack((data, np.c_[xs, ys]))
if i == j:
targets.extend([1] * count)
else:
targets.extend([-1] * count)
return np.c_[data, targets]
1logistic_regression.py 文件源码
项目：Machine-Learning-Algorithms
作者: PacktPublishing
项目源码
文件源码
阅读 22
收藏 0
点赞 0
评论 0
def show_classification_areas(X, Y, lr):
x_min, x_max = X[:, 0].min() - .5, X[:, 0].max() + .5
y_min, y_max = X[:, 1].min() - .5, X[:, 1].max() + .5
xx, yy = np.meshgrid(np.arange(x_min, x_max, 0.02), np.arange(y_min, y_max, 0.02))
Z = lr.predict(np.c_[xx.ravel(), yy.ravel()])
Z = Z.reshape(xx.shape)
plt.figure(1, figsize=(30, 25))
plt.pcolormesh(xx, yy, Z, cmap=plt.cm.Pastel1)
# Plot also the training points
plt.scatter(X[:, 0], X[:, 1], c=np.abs(Y - 1), edgecolors='k', cmap=plt.cm.coolwarm)
plt.xlabel('X')
plt.ylabel('Y')
plt.xlim(xx.min(), xx.max())
plt.ylim(yy.min(), yy.max())
plt.xticks(())
plt.yticks(())
plt.show()
def plot_decision_boundary(X, Y, model):
# X - some data in 2dimensional np.array
x_min, x_max = X[:, 0].min() - 1, X[:, 0].max() + 1
y_min, y_max = X[:, 1].min() - 1, X[:, 1].max() + 1
xx, yy = np.meshgrid(np.arange(x_min, x_max, 0.01),
np.arange(y_min, y_max, 0.01))
# here "model" is your model's prediction (classification) function
Z = model(np.c_[xx.ravel(), yy.ravel()])
# Put the result into a color plot
Z = Z.reshape(xx.shape)
plt.contourf(xx, yy, Z, cmap=plt.cm.Paired)
plt.axis('off')
for i in x:
print i
# Plot also the training points
plt.scatter(X[:, 0], X[:, 1], c=Y, cmap=plt.cm.Paired)
#???????
def get_double_integrator(m=1000, b=50, d=1):
N = 2
sys = LTISystem(
np.c_[[0, 1], [1, -b/m]], # A
np.r_[0, d/m], # B
# np.r_[0, 1], # C
)
def ref_func(*args):
if len(args) == 1:
x = np.zeros(N)
else:
x = args[1]
return np.r_[d/m, 0]-x
ref = SystemFromCallable(ref_func, N, N)
return sys, ref
def get_electromechanical(b=1, R=1, L=1, K=np.pi/5, M=1):
# TODO: determine good reference and/or initial_condition
# TODO: determine good default values for b, R, L, M
N = 3
sys = LTISystem(
np.c_[ # A
[0, 0, 0],
[1, -b/M, -K/L],
[0, K/M, -R/L]
],
np.r_[0, 0, 1/L], # B
# np.r_[1, 0, 0], # C
)
sys.initial_condition = np.ones(N)
def ref_func(*args):
if len(args) == 1:
x = np.zeros(N)
else:
x = args[1]
return np.r_[0, 0, 0]-x
ref = SystemFromCallable(ref_func, N, N)
return sys, ref
def get_cart_pendulum(m=1, M=3, L=0.5, g=9.81, pedant=False):
N = 4
sys = LTISystem(
np.c_[ # A
[0, 0, 0, 0],
[1, 0, 0, 0],
[0, m*g/M, 0, (-1)**(pedant)*(m+M)*g/(M*L)],
[0, 0, 1, 0]
],
np.r_[0, 1/M, 0, 1/(M*L)], # B
# np.r_[1, 0, 1, 0] # C
)
sys.initial_condition = np.r_[0, 0, np.pi/3, 0]
def ref_func(*args):
if len(args) == 1:
x = np.zeros(N)
else:
x = args[1]
return np.zeros(N)-x
ref = SystemFromCallable(ref_func, N, N)
return sys, ref
def update_equation_function(self, *args):
event_var = self.event_variable_equation_function(*args)
if self.condition_idx is None:
self.condition_idx = np.where(np.all(np.r_[
np.c_[[[True]], event_var >= self.event_bounds],
np.c_[event_var <= self.event_bounds, [[True]]]
], axis=0))[0][0]
else:
sq_dist = (event_var - self.event_bounds)**2
crossed_root_idx = np.where(sq_dist == np.min(sq_dist))[1][0]
if crossed_root_idx == self.condition_idx:
self.condition_idx += 1
elif crossed_root_idx == self.condition_idx-1:
self.condition_idx -= 1
else:
warnings.warn("SwitchedSystem did not cross a neighboring " +
"boundary. This may indicate an integration " +
"error. Continuing without updating " +
"condition_idx", UserWarning)
return self.state_update_equation_function(*args)
def plane2xyz(center, ij, plane):
"""
converts image pixel indices to xyz on the PLANE.
center : 2-tuple
ij : nx2 int array
plane : 4-tuple
return nx3 array.
"""
ij = np.atleast_2d(ij)
n = ij.shape[0]
ij = ij.astype('float')
xy_ray = (ij-center[None,:]) / DepthCamera.f
z = -plane[2]/(xy_ray.dot(plane[:2])+plane[3])
xyz = np.c_[xy_ray, np.ones(n)] * z[:,None]
return xyz
