用python画动态三维轨迹_如何用matplotlib绘制python中的3D密度图

感谢mwaskon – 建议mayavi图书馆.

我在mayavi中重新创建了密度散点图,如下所示:

import numpy as np

from scipy import stats

from mayavi import mlab

mu, sigma = 0, 0.1

x = 10*np.random.normal(mu, sigma, 5000)

y = 10*np.random.normal(mu, sigma, 5000)

z = 10*np.random.normal(mu, sigma, 5000)

xyz = np.vstack([x,y,z])

kde = stats.gaussian_kde(xyz)

density = kde(xyz)

# Plot scatter with mayavi

figure = mlab.figure('DensityPlot')

pts = mlab.points3d(x, y, z, density, scale_mode='none', scale_factor=0.07)

mlab.axes()

mlab.show()

将scale_mode设置为“none”可防止与密度向量成比例缩放字形.另外对于大型数据集,我禁用场景渲染,并使用掩码减少点数.

# Plot scatter with mayavi

figure = mlab.figure('DensityPlot')

figure.scene.disable_render = True

pts = mlab.points3d(x, y, z, density, scale_mode='none', scale_factor=0.07)

mask = pts.glyph.mask_points

mask.maximum_number_of_points = x.size

mask.on_ratio = 1

pts.glyph.mask_input_points = True

figure.scene.disable_render = False

mlab.axes()

mlab.show()

接下来,评估网格上的高斯kde:

import numpy as np

from scipy import stats

from mayavi import mlab

mu, sigma = 0, 0.1

x = 10*np.random.normal(mu, sigma, 5000)

y = 10*np.random.normal(mu, sigma, 5000)

z = 10*np.random.normal(mu, sigma, 5000)

xyz = np.vstack([x,y,z])

kde = stats.gaussian_kde(xyz)

# Evaluate kde on a grid

xmin, ymin, zmin = x.min(), y.min(), z.min()

xmax, ymax, zmax = x.max(), y.max(), z.max()

xi, yi, zi = np.mgrid[xmin:xmax:30j, ymin:ymax:30j, zmin:zmax:30j]

coords = np.vstack([item.ravel() for item in [xi, yi, zi]])

density = kde(coords).reshape(xi.shape)

# Plot scatter with mayavi

figure = mlab.figure('DensityPlot')

grid = mlab.pipeline.scalar_field(xi, yi, zi, density)

min = density.min()

max=density.max()

mlab.pipeline.volume(grid, vmin=min, vmax=min + .5*(max-min))

mlab.axes()

mlab.show()

作为最后的改进,我通过并行调用kde函数加快了密度密度函数的评估.

import numpy as np

from scipy import stats

from mayavi import mlab

import multiprocessing

def calc_kde(data):

return kde(data.T)

mu, sigma = 0, 0.1

x = 10*np.random.normal(mu, sigma, 5000)

y = 10*np.random.normal(mu, sigma, 5000)

z = 10*np.random.normal(mu, sigma, 5000)

xyz = np.vstack([x,y,z])

kde = stats.gaussian_kde(xyz)

# Evaluate kde on a grid

xmin, ymin, zmin = x.min(), y.min(), z.min()

xmax, ymax, zmax = x.max(), y.max(), z.max()

xi, yi, zi = np.mgrid[xmin:xmax:30j, ymin:ymax:30j, zmin:zmax:30j]

coords = np.vstack([item.ravel() for item in [xi, yi, zi]])

# Multiprocessing

cores = multiprocessing.cpu_count()

pool = multiprocessing.Pool(processes=cores)

results = pool.map(calc_kde, np.array_split(coords.T, 2))

density = np.concatenate(results).reshape(xi.shape)

# Plot scatter with mayavi

figure = mlab.figure('DensityPlot')

grid = mlab.pipeline.scalar_field(xi, yi, zi, density)

min = density.min()

max=density.max()

mlab.pipeline.volume(grid, vmin=min, vmax=min + .5*(max-min))

mlab.axes()

mlab.show()

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值