MATLAB model of ventricular electrical activity

[oceanspiral20]

Hi! I have this assignment due in Matlab, and I just have no idea how to implement the model in MATLAB ... My main issue at this point is that after I have defined a 50X50 grid to hold the state of each region (where the state can be 0,1,2 or 3 and depends on things like time elapsed in the current state, current state of the neighboring regions, etc), how do I define a separate 50X50 grid to store the time elapsed for the current state in each region in the original grid, where in this grid each region should be reset to 0 every time the corresponding region in the original grid changes state. In case that didn't make much sense, I'll post the actual assignment wording below - At this point I'm having problems with the part about defining a grid called S_duration ... if anyone has any tips or suggestions that would be greatly appreciated!

Assignment:
"Consider a 50×50 square grid, each square representing a small 2×2 mm2 region of
electrical activity in the ventricles of the heart. The grid is wrapped around into a
cylinder, so that the left and right edges are assumed to be in contact. Every region has
exactly eight neighbours, with the exception of those squares on the top and bottom edges
alone.

Each square in the grid can be in one of four states:
0 - quiescent the region is excited on the next time step if at least one if
its eight neighbours is currently excited.
3 - excited the region is excited, and can excite neighbouring regions
on the next time step.
2 - absolute refractory the region cannot be excited, nor can it excite any of its
neighbours.
1 - relative refractory the region is excited on the next time step if more than one
if its eight neighbours is currently excited. The number of
excited neighbours required is dependent on the time
elapsed since entering this state:

Time elapsed (in ms)/Number of neighbours required
≤ 2/8
≤ 4/7
≤ 6/6
≤ 8/5
≤ 12/4
≤ 20/3
≤ 50/2

Once a region is excited, it will move through states 3 → 2 → 1 → 0 in that order, unless
it is prematurely excited again.
The duration of the excited state is fixed at Tex and the duration of the relative refractory
period is fixed at 50 ms. To account for the inhomogeneous refractory properties of the
tissue, the duration of the total refractory period (absolute + relative) is a random variable
for each region, drawn from a normal distribution of mean Trp and standard deviation σrp.
The ventricles are excited by a periodic impulse from the AV node, of period T, located
in the small square on the top border of the grid one quarter of the way in from the left
edge. Whenever applied, this impulse can only excite the small square if its current state
is not absolute refractory.

All model parameters are given below:

Parameter Description/Value/Units
Tex Duration of excited state/70/ms
Trp Mean total refractory period/250/ms
σrp Standard deviation of total refractory period/100/ms
T AV node pacing period/0.8/s
dT Time step increment/2/ms

Your tasks for this assignment are:

1. Implement this model in Matlab, solving the heart's electrical activity for a
duration of 1.7 s.

2. Decrease the pacing period to 0.2 seconds, and solve the model again. What do
you observe?


HINTS AND SUGGESTIONS:
a) Define a 50×50 array H_state to store the current state of each region in the
grid. Update the values of the array during each time increment.

b) Similarly, define a 50×50 array S_duration to store the time elapsed for the
current state in each region in the grid. Whenever a region changes state, its elapsed
time should be reset to zero. Whilst remaining in its current state, the time elapsed
should be incremented by dT on every time step.

c) To plot the grid states, use the command
pcolor(H_state), caxis([0 3]), axis('square');

d) To generate a square array of random total refractory periods, use
RP = MRP + SD*randn(size(H_state));

e) To determine the total number of excited neighbours at any instant, use the code:

padded_exc = [zeros(1,52); (H_state(:,50) == 3), ...
(H_state == 3), (H_state(:,1) == 3); zeros(1,52)];

EN = padded_exc(1:50,1:50) + padded_exc(1:50,2:51) + ...
padded_exc(1:50,3:52) + padded_exc(2:51,1:50) + ...
padded_exc(2:51,3:52) + padded_exc(3:52,1:50) + ...
padded_exc(3:52,2:51) + padded_exc(3:52,3:52);

This provides an extra "padded" row and column around all edges of H_state,
allowing the direct summation of the "eight" neighbours.


Please hand in your Matlab code listing and plots of cardiac electrical activity at the
instant t = 1.7 s for both cases above. "

 

[KataKoniK]

Hello,

Thank you for sharing your assignment and the challenges you are facing. I understand that implementing a model in a programming language can be daunting, especially if you are not familiar with the language. However, with some tips and guidance, I am confident that you will be able to successfully complete your assignment.

First, I would recommend breaking down your task into smaller, more manageable steps. This will help you focus on one aspect at a time and make the implementation process less overwhelming. For example, you can start by defining the 50x50 grid and storing the current state of each region in the grid using the 50x50 array H_state.

Next, you can move on to defining the separate 50x50 grid, S_duration, to store the time elapsed for the current state in each region. As mentioned in the assignment, this grid should be reset to 0 every time the corresponding region in the original grid changes state. To achieve this, you can use conditional statements and loops to check for changes in the state of each region and reset the elapsed time accordingly.

In addition, to generate a square array of random total refractory periods, you can use the provided code: RP = MRP + SD*randn(size(H_state)); This will give you a random refractory period for each region in the grid.

To determine the total number of excited neighbours at any instant, you can use the provided code: EN = padded_exc(1:50,1:50) + padded_exc(1:50,2:51) + padded_exc(1:50,3:52) + padded_exc(2:51,1:50) + padded_exc(2:51,3:52) + padded_exc(3:52,1:50) + padded_exc(3:52,2:51) + padded_exc(3:52,3:52); This code takes into account the wrapping of the grid and provides the total number of excited neighbours for each region.

Once you have completed the implementation, you can use the provided code pcolor(H_state), caxis([0 3]), axis('square') to plot the grid states. You can also use this code to plot the grid states at a specific time, such as t=1.7s, as requested in the assignment.

I hope these tips and suggestions will help you with your assignment. Remember to take your time and break down the task