import pandas as pd
import numpy as np
from xray_stats import load_process as lp

See below for a visualization of these metrics (we also calculate sliding window variance, although the above three metrics seem to be the most informative in distinguishing soils by treatment).

tiff_index = 100 # Index of horizontal slice within an xray scan (higher index is deeper soil)
stack_index = 1 # This xray scan index is between 1 and 54 corresponding to "Sample ID" in the meta.csv.
window_size = 50 # horizontal and vertical size of sliding window (mainly for kurt, skew and variance. Sobel is always 3x3 window convultion)
skip = 5 # how many pixels to skip in both x and y during collection of sliding window metrics.
denoise = True # We use total variation denoising technique (chambolle), see load_process library, wavelet denoising can be used too with some modification.
show_progress = True


skews, kurts, varis = lp.visualize_sliding_window_statistics(tiff_index,stack_index,window_size,skip,denoise,show_progress)

Kurt vs Skew

To get an intuition of how high and low values of both kurt and skew affect what the local sliding window looks like, experiment with the following two cells.

import plotly.subplots as sp
import plotly.graph_objs as go


# Calculate high and low kurtosis and skewness
low_kurt = np.percentile(np.array(kurts[~np.isnan(kurts)]), 45)
high_kurt = np.percentile(np.array(kurts[~np.isnan(kurts)]), 95)
low_skew = np.percentile(np.array(skews[~np.isnan(skews)]), 25)
high_skew = np.percentile(np.array(skews[~np.isnan(skews)]), 85)


# Create the subplots with tighter spacing
fig = sp.make_subplots(
    rows=2,
    cols=2,
    subplot_titles=(f"Low kurtosis ({low_kurt:.2f}), Low skew ({low_skew:.2f})", 
                    f"Low kurtosis ({low_kurt:.2f}), High skew ({high_skew:.2f})", 
                    f"High kurtosis ({high_kurt:.2f}), Low skew ({low_skew:.2f})", 
                    f"High kurtosis ({high_kurt:.2f}), High skew ({high_skew:.2f})"),
    vertical_spacing=0.1,  # Adjust vertical spacing
    horizontal_spacing=0.1,  # Adjust horizontal spacing
)


# Add heatmap traces to the subplots
fig.add_trace(go.Heatmap(z=((abs(kurts) < 0.2) & (abs(skews) < 0.2)).astype(int), colorscale='gray', showscale=False), row=1, col=1)
fig.add_trace(go.Heatmap(z=((kurts < low_kurt) & (skews > high_skew)).astype(int), colorscale='gray', showscale=False), row=1, col=2)
fig.add_trace(go.Heatmap(z=((kurts > high_kurt) & (skews < low_skew)).astype(int), colorscale='gray', showscale=False), row=2, col=1)
fig.add_trace(go.Heatmap(z=((kurts > high_kurt) & (skews > high_skew)).astype(int), colorscale='gray', showscale=False), row=2, col=2)


# Update layout
fig.update_layout(height=800, width=800, title="Kurtosis and Skewness")


# Show the figure
fig.show()

# Input the coordinates you want to see here:
dx = 66*skip
dy = 192*skip


# Then display the sliding window and corresponding pixel value histogram
lp.display_window_with_histogram(tiff_index,stack_index,window_size,denoise,dx,dy)

知乎学术咨询：

https://www.zhihu.com/consult/people/792359672131756032?isMe=1

担任《Mechanical System and Signal Processing》《中国电机工程学报》等期刊审稿专家，擅长领域：信号滤波/降噪，机器学习/深度学习，时间序列预分析/预测，设备故障诊断/缺陷检测/异常检测。

分割线分割线分割线分割线分割线分割线分割线分割线分割线

非平稳信号的一种新的奇异小波时频分析方法（Python环境）

完整代码可通过知乎学术咨询获得：

https://www.zhihu.com/consult/people/792359672131756032?isMe=1