Here is a sample code that reads the data from a text file, filters it, displays it, computes m-wave and h-reflex amplitudes, and plots them:
% Read data from text file
data = importdata('data.txt');
% Set variables
Fs = 1000; % Sampling frequency
t = (0:length(data)-1)/Fs; % Time vector
% Filter data
[b,a] = butter(4, [20, 500]/(Fs/2), 'bandpass');
data_filt = filtfilt(b, a, data);
% Display raw and filtered data
figure;
subplot(2,1,1); plot(t, data); xlabel('Time (s)'); ylabel('Amplitude'); title('Raw data');
subplot(2,1,2); plot(t, data_filt); xlabel('Time (s)'); ylabel('Amplitude'); title('Filtered data');
% Compute m-wave and h-reflex
pre_stim_time = 0.05; % Time before stimulus onset
post_stim_time = 0.2; % Time after stimulus onset
pre_stim_samples = round(pre_stim_time * Fs);
post_stim_samples = round(post_stim_time * Fs);
stim_time = 0.01; % Stimulus duration
stim_samples = round(stim_time * Fs);
threshold = 0.1; % Threshold for peak detection
m_wave_amp = zeros(size(data, 2), 1);
h_reflex_amp = zeros(size(data, 2), 1);
for i = 1:size(data, 2)
% Find stimulus onset
onset = find(data_filt(pre_stim_samples+1:end, i) > threshold, 1) + pre_stim_samples;
if isempty(onset)
continue;
end
% Find m-wave peak
[m_wave_amp(i), m_wave_peak] = max(data_filt(onset:onset+stim_samples, i));
m_wave_peak = m_wave_peak + onset - 1;
% Find h-reflex peak
[h_reflex_amp(i), h_reflex_peak] = max(data_filt(onset+stim_samples:onset+stim_samples+post_stim_samples, i));
h_reflex_peak = h_reflex_peak + onset + stim_samples - 1;
% Plot m-wave and h-reflex
figure;
plot(t(onset-pre_stim_samples:onset+stim_samples+post_stim_samples), data_filt(onset-pre_stim_samples:onset+stim_samples+post_stim_samples, i));
hold on;
plot(t(m_wave_peak), m_wave_amp(i), 'rv', 'MarkerSize', 10, 'LineWidth', 2);
plot(t(h_reflex_peak), h_reflex_amp(i), 'r^', 'MarkerSize', 10, 'LineWidth', 2);
xlabel('Time (s)'); ylabel('Amplitude'); title(sprintf('M-wave = %.2f mV, H-reflex = %.2f mV', m_wave_amp(i), h_reflex_amp(i)));
legend('Filtered data', 'M-wave', 'H-reflex');
end
Before running the code, you need to replace 'data.txt'
with the filename of your text file. The code assumes that the file contains one column of data for each channel. If your file has a different format, you may need to modify the code accordingly.
The code filters the data using a 4th-order Butterworth bandpass filter with cutoff frequencies of 20 Hz and 500 Hz. You can adjust these values to match your data.
The code computes m-wave and h-reflex amplitudes for each channel of data, and plots them along with the filtered data. You can adjust the parameters pre_stim_time
, post_stim_time
, stim_time
, and threshold
to match your stimulus and data.