Hemodynamic Traveling Waves in Human Visual Cortex

Functional MRI (fMRI) experiments rely on precise characterization of the blood oxygen level dependent (BOLD) signal. As the spatial resolution of fMRI reaches the sub-millimeter range, the need for quantitative modelling of spatiotemporal properties of this hemodynamic signal has become pressing. Here, we find that a detailed physiologically-based model of spatiotemporal BOLD responses predicts traveling waves with velocities and spatial ranges in empirically observable ranges. Two measurable parameters, related to physiology, characterize these waves: wave velocity and damping rate. To test these predictions, high-resolution fMRI data are acquired from subjects viewing discrete visual stimuli. Predictions and experiment show strong agreement, in particular confirming BOLD waves propagating for at least 5–10 mm across the cortical surface at speeds of 2–12 mm s-1. These observations enable fundamentally new approaches to fMRI analysis, crucial for fMRI data acquired at high spatial resolution.


SupplementaryMaterial TextS1
Spatiotemporalhemodynamics Themodelequationsdescribedin [1]detailthegeneralhemodynamicresponse.Inthepresent study, the linear spatiotemporal hemodynamic response function (stHRF) is derived from this general model. Further mathematical analysis will be published in a subsequent paper. This document references equations from the Methods section of the main text by number, while equationsfromthesupportinginformationareprefixedS.
The source and sink each correspond to boundary conditions on the fluid velocity terms that enterthemassandmomentumconservationequations.Byconsideringthefluxofbloodvelocity to the mass inflow entering the system, the condition for the inflow fluid velocity v F can be writtenas . (S2) whichisthefinalterminEq.5ofthemaintext.

ModelLinearization
As long as the neural activity signal is sufficiently small, the hemodynamic response can be estimated through linear analysis. Under this assumption one can analyze linear perturbations from the steady state. This is achieved mathematically by writing each variable θ [i.e., either F,Q,P,v,z,orξ]asthesumofitssteadyvalueθ 0 anditslinearperturbationθ 1 ,i.e. .

(S5)
In this system, the steady state is determined by setting all spatial and temporal derivatives to zeroandsolvingtheensuingalgebraicequations.Afurtherassumptionisthatthesteadystateon averageisspatiallyuniformandthatthemeanbloodfluidvelocityaveragesto0.
With these assumptions, the system dynamics can now be represented by four evolution andtheBOLDsignalequationis .
Togetherthesetransferfunctionsforindividualprocessesyieldtheoveralltransferfunctionvia theBOLDsignalEq.S9: . (S16) UsingS13-S15,thetransferfunctionS16fortheBOLDsignalequationcanberewrittenas . (S18) where δ is the Dirac delta function and k x and k y' are the spatial frequencies in the x and y' directions.ThereforetheonedimensionalspatiotemporalHRFtoa1Dneuralstimulusisgiven wherethedeltafunctioninEq.S20hasbeenusedtoevaluatethey'transform.Eq.S21yieldsthe response to a line stimulus on the cortex, such as that evoked by an isoeccentric curve in the visualfield. Torepresenttheexperiment,thefollowingformforzwasused:

Modelvariablesandparameters
The model contains physiological variables and parameters. a,c TableS1:Themodelvariablesandparameters.Ateachrow,thefirstcolumndetailsthequantity, the second column show the symbols that describe these quantities in the model. The third columndetailsthenominalrange(ifitexistsintheliterature),withitsappropriateunitsinthe fourth column and its source in the 5 th column. The last column details notes on each variable/parameter:adonothaveadirectanalogueintheballoonmodel,bisaparameterused inpreviousballoonmodels,carethesetofindependentvariablesneededtocalculatethestHRF.

Rangesfordampingandpropagationvelocity
The theory shows that calculation of spatiotemporal properties of the response requires two parametersinadditiontothosepresentintheballoonmodel,forexample.Aprioriconstraintson these parameters, the propagation velocity ν β and temporal damping rate Γ, can be made by comparisonswithpreviousexperimentalwork.Thismeansthat,althoughpriorworkdoesnot giveprecisevaluesforthesequantities,itdoesyieldtheirapproximatevaluesandpreventsthem frombeingtreatedasfreeparameters.
The temporal damping rate Γ has two components, the average viscous damping and the contributionfromlossduetooutflow,with .
• Red Blood Cell velocity, is in the order of 10 mm s-1 [3], which is attributed to the effectivevelocityv.

[S29]
where n is the degree of the polynomial, N is the total number of points fitted, and σ n is the standarddeviationoftheresidualerrorofthefitfordegreen.