求推荐R中绘制图像堆栈的高效优雅实现方法

我尝试可视化3D数组切片（即图像栈），但发现没有R包提供能轻松完成该任务的内置函数。我基于模拟数据写了一份自定义实现的示例代码，但处理大型CT扫描数据时机器负载极高，想找更简便、计算成本更低的替代方案？

示例代码 #

library(leaflet) # 颜色调色板
library(rgl) # 3D可视化

# 示例数据
slides = array(data = 0,
               dim = c(200,200,15))

# 添加强度值
for(i in 1:dim(slides)[3]) slides[20:180,20:180,i] = i*10

# x/y/z方向的体素尺寸
dims = c(1,1,20) # 调整z轴尺寸用于可视化

# 生成所有可能的x/y/z组合
xyz = expand.grid(1:dim(slides)[1],1:dim(slides)[2],1:dim(slides)[3]) 

# 获取每个数据点的强度值
value = apply(xyz, 1, function(x) slides[x[1],x[2],x[3]]) 
# 反转值
value = max(value) - value 

# 应用体素尺寸到x/y/z
xyz = t(t(xyz) * dims) 

# 将强度值转换为颜色
colb = colorNumeric(palette = "Greys",
                   domain = c(0,max(value)))
col = colb(value)

# 绘制图像栈
plot3d(xyz, col = col, aspect = FALSE, axes = FALSE, xlab = "", ylab = "", zlab = "")  

当前可视化结果 #

高效替代方案 #

1. 用rgl内置体素函数优化 #

你当前代码是把每个体素拆成单独点绘制，数据量大时效率极低。rgl自带的voxels()函数可以直接处理3D数组，无需手动展开坐标：

library(rgl)
library(leaflet)

colb = colorNumeric(palette = "Greys", domain = c(0, max(slides)))

# 直接渲染体素块，匹配尺寸缩放
voxels(slides, 
       color = colb(max(slides) - slides), 
       aspect = dims,
       axes = FALSE, xlab = "", ylab = "", zlab = "")

这种方式直接渲染体素网格，内存占用和计算量会大幅降低。

2. 用rayshader做硬件加速渲染 #

rayshader针对体数据做了性能优化，支持硬件加速，处理大型CT数据时速度优势明显：

library(rayshader)
library(leaflet)

# 调整轴顺序适配rayshader要求
slides_shader = aperm(slides, c(2,1,3))
colb = colorNumeric(palette = "Greys", domain = c(0, max(slides)))
color_map = colb(max(slides) - slides_shader)

# 渲染3D场景，匹配z轴缩放比例
plot_3d(slides_shader, 
        zscale = dims[3]/dims[1],
        solidcolor = color_map,
        axes = FALSE, xlab = "", ylab = "", zlab = "")

3. 用ITKR处理医学CT数据 #

如果是专业医学CT数据，ITKR基于ITK库，专门针对医学影像优化，支持直接加载DICOM序列，内置高效3D渲染：

library(ITKR)

# 加载CT文件夹中的DICOM序列
ct_image = readDicomSeries("你的CT数据文件夹路径")
# 一键可视化3D体数据
view3D(ct_image)

4. 先做阈值过滤减少渲染量 #

如果只需要看特定强度范围的体素，先过滤掉背景无效数据：

library(rgl)
library(leaflet)

# 只保留强度大于50的体素
mask = slides > 50
slides_filtered = slides
slides_filtered[!mask] = NA

colb = colorNumeric(palette = "Greys", domain = c(0, max(slides)))
voxels(slides_filtered, 
       color = colb(max(slides) - slides_filtered),
       aspect = dims,
       axes = FALSE)

这种方式能避免渲染大量无意义的背景体素，直接降低机器负载。

内容的提问来源于stack exchange，提问作者mri