A Mathematica Companion to: 
Systems Factorial Technology Explained to Humans.

(c) Harding, B, Goulet, M.-A., Jolin, S., Tremblay, C., Villeneuve, S., Durand, G.

The purpose of this Mathematica notebook is to show the differences between the five possible processing architectures using the exact same simulation structure. The only change between simulations is the aggregation function of the two sub-processes’ completion times. 

It is therefore intended to be used as a compliment to “Systems Factorial Technology Explained to Humans” to help visualize the various architectures and their characteristics. We hope that this Mathematica notebook will help users understand the mechanics behind the SFT methodology and see why this advancement in mathematica psychology is worth adding to one’s repertoire of analyses.

For the following architectures, we have used the following aggregation functions:
- Serial self-terminating	ChooseOne	i.e., one sub-processes’ RT is selected at random
- Serial exhaustive		Plus		i.e., take the sum of both sub-processes RT
- Parallel self-terminating	Min		i.e., take the fastest RT between both sub-processes
- Parallel exhaustive		Max		i.e., take the slowest RT between both sub-processes
- Coactive			Coactive	i.e., take many incoming times for each process and extract the (k^th) one

General definitions			𝒟, μ, σ, n

Here we define the mean (μ) and standard deviation (σ) for the simulated sub-processes’ response times (RT). All distributions follow a normal distribution:

𝒟[μ_,σ_]:=NormalDistribution[μ,σ];

(Subscript[μ, h1])=20; (*Average time for sub-process "A" in High condition, or, Alice's average time when she is rested*)
(Subscript[μ, l1])=36; (*Average time for sub-process "A" in Low condition, or, Alice's average time when she is tired*)
(Subscript[μ, h2])=23; (*Average time for sub-process "B" in High condition, or, Bob's average time when he is rested*)
(Subscript[μ, l2])=40; (*Average time for sub-process "B" in Low condition, or, Bob's average time when he is tired*)σ=5;

Here we set he number of simulated RT. For clean results, such as the ones presented here, we suggest using a large n. For more realistic results, we encouage using a much smaler n (such as n = 500):

n=50000;

1- Serial Self-Terminating		ChooseOne

Define aggregate function

We must clear all parameter values before use. This ensures that each architecture simulated here is independant and not contaminated with any of the other values.

ClearAll[time1,time2,RThh,RThl,RTlh,RTll,RTh1,RTh2,Hab,Ha,Hb,Kab,Ka,Kb];

This aggregation function randomly selects which sub-processes’ RT will be considered.

ChooseOne[time1_,time2_,p_]:=If[Random[]<p,time1,time2];

Make Distributions

Here we create the RT distributions for each of the four experimental conditions. In these cells “time1” and ”time2” represent the individual completion times of both sub-processes.

time1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],n];
RThh=MapThread[ChooseOne[#1,#2,1/2]&,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, l2]),σ],n];
RThl=MapThread[ChooseOne[#1,#2,1/2]&,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, l1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],n];
RTlh=MapThread[ChooseOne[#1,#2,1/2]&,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, l1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, l2]),σ],n];
RTll=MapThread[ChooseOne[#1,#2,1/2]&,{time1,time2}];

Here we plot all four experimental conditions’ RT on the same plot

SmoothHistogram[{RThh,RThl,RTlh,RTll},Automatic,"PDF",AxesLabel->{"RT","Frequency"}]

Each experimental condition’s RT can be individually plotted. For example, if one would want to plot the HH condition:

Histogram[RThh,Automatic,"PDF"]

MIC with statistical test

MIC

Here we measure the MIC with:

MIC=(Mean[RThh]+Mean[RTll])-(Mean[RThl]+Mean[RTlh])

Statistical tests

And here we test the null hyposis that MIC is not equal to 0.

SE[data_]:=((StandardDeviation[data])/(Sqrt[Length[data]]))

SEmic=(1/2)(Sqrt[((SE[RTll])^2)+((SE[RTlh])^2)+((SE[RThl])^2)+((SE[RThh])^2)]);

pmic=2(1-CDF[NormalDistribution[],Abs[MIC/SEmic]])

Survival Functions

Each of the experimental conditions can be transformed into survivor functions and plotted together.

hi=Max[Flatten[{{RThh,RThl,RTlh,RTll}}]];
bins={0,hi,1};

shh=HistogramList[RThh,bins,"SF"][[2]];
shl=HistogramList[RThl,bins,"SF"][[2]];
slh=HistogramList[RTlh,bins,"SF"][[2]];
sll=HistogramList[RTll,bins,"SF"][[2]];

ListPlot[{shh,shl,slh,sll},Joined->True,AxesLabel->{"RTs","Frequency"}]

SIC Curve with statistical tests

SIC with centerline

Here we measure the centerline proposed by Harding, Leblanc, Goulet, & Cousineau (submitted). It is the median of all experimental conditions pooled together.

centerline=Median[Flatten[{{RThh,RThl,RTlh,RTll}}]];

Now we can measure and plot the SIC (along with the SIC centerline) with:

SIC=(shh+sll)-(shl+slh);

ListPlot[SIC,Joined->True,PlotRange->{All,{-1,1}},AxesLabel->{"t ","SIC(t)"},GridLines->{{{centerline,{Red,Dashed}}}}]

Kolmogorov-Smirnov test

Here we use an adapted Kolmogorov-Smirnov test (Houpt & Townsend, 2010) to see whether the SIC differs from 0. As we already have the SIC, which is shown above, all we need is the Harmonic Mean.

hMean=HarmonicMean[{Length[RThh],Length[RTlh],Length[RThl],Length[RTll]}];

Here we test the positive portion of the curve

(SuperPlus[D])=Max[SIC]
(SuperPlus[P])=Exp[-2(hMean/4)((((SuperPlus[D])))^2)]

Here we test the negative portion of the curve

(SuperMinus[D])=Min[SIC]
(SuperMinus[P])=Exp[-2(hMean/4)((((SuperMinus[D])))^2)]

Capacity

Here we test the capacity of our cognitive process. First, we must gather RT when a sub-process is working alone. In these cells, “RTh1”, “RTh2” are the individual RT for both sub-processes working alone in their H condition and “RThh” is the systems’s RT when both sub-processes are operating together.

RTh1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],n];
RTh2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],n];
RThh=MapThread[ChooseOne[#1,#2,1/2]&,{RTh1,RTh2}];

Finally, using the cumulative hazard function of each of the three working conditions (RTh1, RTh2, RThh) we are able to calculate and plot the capacity curve using CapOR as introduced in Townsend & Nozawa (1995) and CapOR1, a variant of CapOR that is shown in Houpt, Blaha, McIntire, Havig, & Townsend (2014). We use CapOR as we are working in a self-terminating environment.

hi=Max[Flatten[{{RThh,RTh1,RTh2}}]];
bins={0,hi,1};

Hab=HistogramList[RThh,bins,"CHF"][[2]];
Ha=HistogramList[RTh1,bins,"CHF"][[2]];
Hb=HistogramList[RTh2,bins,"CHF"][[2]];

Off[Power::infy];
Off[Infinity::indet];

CapOR=(Hab/((Ha+Hb)));
CapOR1=Hab-(Ha+Hb);


ListPlot[CapOR,Joined->True,PlotRange->{All,All},AxesLabel->{"t ","C_(OR)"}]
ListPlot[CapOR1,Joined->True,PlotRange->{All,All},AxesLabel->{"t ","C_(OR)"}]

The “Off” command above ensures that the calculations keep going regardless of any (1/0)calculations that may be present.

2- Serial Exhaustive			Plus

Define aggregate function

We must clear all parameter values before use. This ensures that each architecture simulated here is independant and not contaminated with any of the other values.

ClearAll[time1,time2,RThh,RThl,RTlh,RTll,RTh1,RTh2,Hab,Ha,Hb,Kab,Ka,Kb];

There is no need to define an aggregate function. RT is dictated in the next section by “Plus”.

Make Distributions

Here we create the RT distributions for each of the four experimental conditions. In these cells “time1” and ”time2” represent the individual completion times of both sub-processes.

time1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],n];
RThh=MapThread[Plus,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, l2]),σ],n];
RThl=MapThread[Plus,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, l1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],n];
RTlh=MapThread[Plus,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, l1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, l2]),σ],n];
RTll=MapThread[Plus,{time1,time2}];

Here we plot all four experimental conditions’ RT on the same plot

SmoothHistogram[{RThh,RThl,RTlh,RTll},Automatic,"PDF",AxesLabel->{"RT","Frequency"}]

MIC with statistical test

MIC

Here we measure the MIC with:

MIC=(Mean[RThh]+Mean[RTll])-(Mean[RThl]+Mean[RTlh])

Statistical tests

And here we test the null hyposis that MIC is not equal to 0.

SE[data_]:=((StandardDeviation[data])/(Sqrt[Length[data]]))

SEmic=(1/2)(Sqrt[((SE[RTll])^2)+((SE[RTlh])^2)+((SE[RThl])^2)+((SE[RThh])^2)]);

pmic=2(1-CDF[NormalDistribution[],Abs[MIC/SEmic]])

Survival Functions

Each of the experimental conditions can be transformed into survivor functions and plotted together.

hi=Max[Flatten[{{RThh,RThl,RTlh,RTll}}]];
bins={0,hi,1};

shh=HistogramList[RThh,bins,"SF"][[2]];
shl=HistogramList[RThl,bins,"SF"][[2]];
slh=HistogramList[RTlh,bins,"SF"][[2]];
sll=HistogramList[RTll,bins,"SF"][[2]];

ListPlot[{shh,shl,slh,sll},Joined->True,AxesLabel->{"RTs","Frequency"}]

SIC Curve with statistical tests

SIC with centerline

Here we measure the centerline proposed by Harding, Leblanc, Goulet, & Cousineau (submitted). It is the median of all experimental conditions pooled together.

centerline=Median[Flatten[{{RThh,RThl,RTlh,RTll}}]];

Now we can measure and plot the SIC (along with the SIC centerline) with:

SIC=(shh+sll)-(shl+slh);

ListPlot[SIC,Joined->True,PlotRange->{All,{-1,1}},AxesLabel->{"t ","SIC(t)"},GridLines->{{{centerline,{Red,Dashed}}}}]

Kolmogorov-Smirnov test

Here we use an adapted Kolmogorov-Smirnov test (Houpt & Townsend, 2010) to see whether the SIC differs from 0. As we already have the SIC, which is shown above, all we need is the Harmonic Mean.

hMean=HarmonicMean[{Length[RThh],Length[RTlh],Length[RThl],Length[RTll]}];

Here we test the positive portion of the curve

(SuperPlus[D])=Max[SIC]
(SuperPlus[P])=Exp[-2(hMean/4)((((SuperPlus[D])))^2)]

Here we test the negative portion of the curve

(SuperMinus[D])=Min[SIC]
(SuperMinus[P])=Exp[-2(hMean/4)((((SuperMinus[D])))^2)]

Capacity

Here we test the capacity of our cognitive process. First, we must gather RT when a sub-process is working alone. In these cells, “RTh1”, “RTh2” are the individual RT for both sub-processes working alone in their H condition and “RThh” is the systems’s RT when both sub-processes are operating together.

RTh1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],n];
RTh2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],n];
RThh=MapThread[Plus,{RTh1,RTh2}];

Finally, using the reverse cumulative hazard function of each of the three working conditions (RTh1, RTh2, RThh) we are able to calculate and plot the capacity curve using CapAND as introduced in Townsend & Nozawa (1995) and CapAND1, a variant of CapAND that is shown in Houpt, Blaha, McIntire, Havig, & Townsend (2014). We use CapAND as we are working in an exhaustive environment.

hi=Max[Flatten[{{RThh,RTh1,RTh2}}]];
bins={0,hi,1};

Kab=HistogramList[RThh,bins,"CDF"][[2]];
Ka=HistogramList[RTh1,bins,"CDF"][[2]];
Kb=HistogramList[RTh2,bins,"CDF"][[2]];
{Ka,Kb,Kab}={Ka,Kb,Kab}//Log;

Off[Power::infy];
Off[Infinity::indet];

CapAND=(((Ka+Kb))/Kab);
CapAND1=Kab-(Ka+Kb);


ListPlot[CapAND,Joined->True,PlotRange->{All,All},AxesLabel->{"t ","C_(AND)"}]
ListPlot[CapAND1,Joined->True,PlotRange->{All,All},AxesLabel->{"t ","C_(AND)"}]

The “Off” command above ensures that the calculations keep going regardless of any (1/0)calculations that may be present.

3- Parallel Self-Terminating		Min

Define aggregate function

We must clear all parameter values before use. This ensures that each architecture simulated here is independant and not contaminated with any of the other values.

ClearAll[time1,time2,RThh,RThl,RTlh,RTll,RTh1,RTh2,Hab,Ha,Hb,Kab,Ka,Kb];

There is no need to define an aggregate function. RT is dictated in the next section by “Min”.

Make Distributions

Here we create the RT distributions for each of the four experimental conditions. In these cells “time1” and ”time2” represent the individual completion times of both sub-processes.

time1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],n];
RThh=MapThread[Min,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, l2]),σ],n];
RThl=MapThread[Min,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, l1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],n];
RTlh=MapThread[Min,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, l1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, l2]),σ],n];
RTll=MapThread[Min,{time1,time2}];

Here we plot all four experimental conditions’ RT on the same plot

SmoothHistogram[{RThh,RThl,RTlh,RTll},Automatic,"PDF",AxesLabel->{"RT","Frequency"}]

MIC with statistical test

MIC

Here we measure the MIC with:

MIC=(Mean[RThh]+Mean[RTll])-(Mean[RThl]+Mean[RTlh])

Statistical tests

And here we test the null hyposis that MIC is not equal to 0.

SE[data_]:=((StandardDeviation[data])/(Sqrt[Length[data]]))

SEmic=(1/2)(Sqrt[((SE[RTll])^2)+((SE[RTlh])^2)+((SE[RThl])^2)+((SE[RThh])^2)]);

pmic=2(1-CDF[NormalDistribution[],Abs[MIC/SEmic]])

Survival Functions

Each of the experimental conditions can be transformed into survivor functions and plotted together.

hi=Max[Flatten[{{RThh,RThl,RTlh,RTll}}]];
bins={0,hi,1};

shh=HistogramList[RThh,bins,"SF"][[2]];
shl=HistogramList[RThl,bins,"SF"][[2]];
slh=HistogramList[RTlh,bins,"SF"][[2]];
sll=HistogramList[RTll,bins,"SF"][[2]];

ListPlot[{shh,shl,slh,sll},Joined->True,AxesLabel->{"RTs","Frequency"}]

SIC Curve with statistical tests

SIC with centerline

Here we measure the centerline proposed by Harding, Leblanc, Goulet, & Cousineau (submitted). It is the median of all experimental conditions pooled together.

centerline=Median[Flatten[{{RThh,RThl,RTlh,RTll}}]];

Now we can measure and plot the SIC (along with the SIC centerline) with:

SIC=(shh+sll)-(shl+slh);

ListPlot[SIC,Joined->True,PlotRange->{All,{-1,1}},AxesLabel->{"t ","SIC(t)"},GridLines->{{{centerline,{Red,Dashed}}}}]

Kolmogorov-Smirnov test

Here we use an adapted Kolmogorov-Smirnov test (Houpt & Townsend, 2010) to see whether the SIC differs from 0. As we already have the SIC, which is shown above, all we need is the Harmonic Mean.

hMean=HarmonicMean[{Length[RThh],Length[RTlh],Length[RThl],Length[RTll]}];

Here we test the positive portion of the curve

(SuperPlus[D])=Max[SIC]
(SuperPlus[P])=Exp[-2(hMean/4)((((SuperPlus[D])))^2)]

Here we test the negative portion of the curve

(SuperMinus[D])=Min[SIC]
(SuperMinus[P])=Exp[-2(hMean/4)((((SuperMinus[D])))^2)]

Capacity

Here we test the capacity of our cognitive process. First, we must gather RT when a sub-process is working alone. In these cells, “RTh1”, “RTh2” are the individual RT for both sub-processes working alone in their H condition and “RThh” is the systems’s RT when both sub-processes are operating together.

RTh1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],n];
RTh2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],n];
RThh=MapThread[Min,{RTh1,RTh2}];

Finally, using the cumulative hazard function of each of the three working conditions (RTh1, RTh2, RThh) we are able to calculate and plot the capacity curve using CapOR as introduced in Townsend & Nozawa (1995) and CapOR1, a variant of CapOR that is shown in Houpt, Blaha, McIntire, Havig, & Townsend (2014). We use CapOR as we are working in a self-terminating environment.

hi=Max[Flatten[{{RThh,RTh1,RTh2}}]];
bins={0,hi,1};

Hab=HistogramList[RThh,bins,"CHF"][[2]];
Ha=HistogramList[RTh1,bins,"CHF"][[2]];
Hb=HistogramList[RTh2,bins,"CHF"][[2]];

Off[Power::infy];
Off[Infinity::indet];

CapOR=(Hab/((Ha+Hb)));
CapOR1=Hab-(Ha+Hb);


ListPlot[CapOR,Joined->True,PlotRange->{All,All},AxesLabel->{"t ","C_(OR)"}]
ListPlot[CapOR1,Joined->True,PlotRange->{All,All},AxesLabel->{"t ","C_(OR)"}]

The “Off” command above ensures that the calculations keep going regardless of any (1/0)calculations that may be present.

4- Parallel Exhaustive			Max

Define aggregate function

There is no need to define the aggregate function. RT is dictated in the next section by “Max”.

We must clear all parameter values before use. This ensures that each architecture simulated here is independant and not contaminated with any of the other values.

ClearAll[time1,time2,RThh,RThl,RTlh,RTll,RTh1,RTh2,Hab,Ha,Hb,Kab,Ka,Kb];

There is no need to define an aggregate function. RT is dictated in the next section by “Max”.

Make Distributions

Here we create the RT distributions for each of the four experimental conditions. In these cells “time1” and ”time2” represent the individual completion times of both sub-processes.

time1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],n];
RThh=MapThread[Max,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, l2]),σ],n];
RThl=MapThread[Max,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, l1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],n];
RTlh=MapThread[Max,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, l1]),σ],n];
time2=RandomVariate[𝒟[(Subscript[μ, l2]),σ],n];
RTll=MapThread[Max,{time1,time2}];

Here we plot all four experimental conditions’ RT on the same plot

SmoothHistogram[{RThh,RThl,RTlh,RTll},Automatic,"PDF",AxesLabel->{"RT","Frequency"}]

MIC with statistical test

MIC

Here we measure the MIC with:

MIC=(Mean[RThh]+Mean[RTll])-(Mean[RThl]+Mean[RTlh])

Statistical tests

And here we test the null hyposis that MIC is not equal to 0.

SE[data_]:=((StandardDeviation[data])/(Sqrt[Length[data]]))

SEmic=(1/2)(Sqrt[((SE[RTll])^2)+((SE[RTlh])^2)+((SE[RThl])^2)+((SE[RThh])^2)]);

pmic=2(1-CDF[NormalDistribution[],Abs[MIC/SEmic]])

Survival Functions

Each of the experimental conditions can be transformed into survivor functions and plotted together.

hi=Max[Flatten[{{RThh,RThl,RTlh,RTll}}]];
bins={0,hi,1};

shh=HistogramList[RThh,bins,"SF"][[2]];
shl=HistogramList[RThl,bins,"SF"][[2]];
slh=HistogramList[RTlh,bins,"SF"][[2]];
sll=HistogramList[RTll,bins,"SF"][[2]];

ListPlot[{shh,shl,slh,sll},Joined->True,AxesLabel->{"RTs","Frequency"}]

SIC Curve with statistical tests

SIC with centerline

Here we measure the centerline proposed by Harding, Leblanc, Goulet, & Cousineau (submitted). It is the median of all experimental conditions pooled together.

centerline=Median[Flatten[{{RThh,RThl,RTlh,RTll}}]];

Now we can measure and plot the SIC (along with the SIC centerline) with:

SIC=(shh+sll)-(shl+slh);

ListPlot[SIC,Joined->True,PlotRange->{All,{-1,1}},AxesLabel->{"t ","SIC(t)"},GridLines->{{{centerline,{Red,Dashed}}}}]

Kolmogorov-Smirnov test

Here we use an adapted Kolmogorov-Smirnov test (Houpt & Townsend, 2010) to see whether the SIC differs from 0. As we already have the SIC, which is shown above, all we need is the Harmonic Mean.

hMean=HarmonicMean[{Length[RThh],Length[RTlh],Length[RThl],Length[RTll]}];

Here we test the positive portion of the curve

(SuperPlus[D])=Max[SIC]
(SuperPlus[P])=Exp[-2(hMean/4)((((SuperPlus[D])))^2)]

Here we test the negative portion of the curve

(SuperMinus[D])=Min[SIC]
(SuperMinus[P])=Exp[-2(hMean/4)((((SuperMinus[D])))^2)]

Capacity

Here we test the capacity of our cognitive process. First, we must gather RT when a sub-process is working alone. In these cells, “RTh1”, “RTh2” are the individual RT for both sub-processes working alone in their H condition and “RThh” is the systems’s RT when both sub-processes are operating together.

RTh1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],n];
RTh2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],n];
RThh=MapThread[Max,{RTh1,RTh2}];

Finally, using the reverse cumulative hazard function of each of the three working conditions (RTh1, RTh2, RThh) we are able to calculate and plot the capacity curve using CapAND as introduced in Townsend & Nozawa (1995) and CapAND1, a variant of CapAND that is shown in Houpt, Blaha, McIntire, Havig, & Townsend (2014). We use CapAND as we are working in an exhaustive environment.

hi=Max[Flatten[{{RThh,RTh1,RTh2}}]];
bins={0,hi,1};

Kab=HistogramList[RThh,bins,"CDF"][[2]];
Ka=HistogramList[RTh1,bins,"CDF"][[2]];
Kb=HistogramList[RTh2,bins,"CDF"][[2]];
{Ka,Kb,Kab}={Ka,Kb,Kab}//Log;

Off[Power::infy];
Off[Infinity::indet];

CapAND=(((Ka+Kb))/Kab);
CapAND1=Kab-(Ka+Kb);


ListPlot[CapAND,Joined->True,PlotRange->{All,All},AxesLabel->{"t ","C_(AND)"}]
ListPlot[CapAND1,Joined->True,PlotRange->{All,All},AxesLabel->{"t ","C_(AND)"}]

The “Off” command above ensures that the calculations keep going regardless of any (1/0)calculations that may be present.

5- Coactive				Coactive

Define aggregate function

We must clear all parameter values before use. This ensures that each architecture simulated here is independant and not contaminated with any of the other values.

ClearAll[time1,time2,RThh,RThl,RTlh,RTll,RTh1,RTh2,Hab,Ha,Hb,Kab,Ka,Kb];

Here we define the aggreagation function for a coactive process. In this function, the times for both sub-processes are pooled. Coactive is therefore a collaboration of both sub-processes.

Coactive[time1_,time2_,k_]:=Sort[Join[Accumulate[time1],Accumulate[time2]]][[k]];

Here we define the two properties of a coactive process: the threshold, “k”, and redundancy, ”redun”. Redundancy ensures that there are enough inputs at any given time that are ready to be sampled. When k=10 channels are needed, redundacy ensures that 30 are present

k=10;

Make Distributions

Here we create the RT distributions for each of the four experimental conditions. In these cells “time1” and ”time2” represent the individual completion times of both sub-processes.

time1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],{n,k}];
time2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],{n,k}];
RThh=MapThread[Coactive[#1,#2,k]&,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],{n,k}];
time2=RandomVariate[𝒟[(Subscript[μ, l2]),σ],{n,k}];
RThl=MapThread[Coactive[#1,#2,k]&,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, l1]),σ],{n,k}];
time2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],{n,k}];
RTlh=MapThread[Coactive[#1,#2,k]&,{time1,time2}];

time1=RandomVariate[𝒟[(Subscript[μ, l1]),σ],{n,k}];
time2=RandomVariate[𝒟[(Subscript[μ, l2]),σ],{n,k}];
RTll=MapThread[Coactive[#1,#2,k]&,{time1,time2}];

Here we plot all four experimental conditions’ RT on the same plot

SmoothHistogram[{RThh,RThl,RTlh,RTll},Automatic,"PDF",AxesLabel->{"RT","Frequency"}]

MIC with statistical test

MIC

Here we measure the MIC with:

MIC=(Mean[RThh]+Mean[RTll])-(Mean[RThl]+Mean[RTlh])

Statistical tests

And here we test the null hyposis that MIC is not equal to 0.

SE[data_]:=((StandardDeviation[data])/(Sqrt[Length[data]]))

SEmic=(1/2)(Sqrt[((SE[RTll])^2)+((SE[RTlh])^2)+((SE[RThl])^2)+((SE[RThh])^2)]);

pmic=2(1-CDF[NormalDistribution[],Abs[MIC/SEmic]])

Survival Functions

Each of the experimental conditions can be transformed into survivor functions and plotted together.

hi=Max[Flatten[{{RThh,RThl,RTlh,RTll}}]];
bins={0,hi,1};

shh=HistogramList[RThh,bins,"SF"][[2]];
shl=HistogramList[RThl,bins,"SF"][[2]];
slh=HistogramList[RTlh,bins,"SF"][[2]];
sll=HistogramList[RTll,bins,"SF"][[2]];

ListPlot[{shh,shl,slh,sll},Joined->True,AxesLabel->{"RTs","Frequency"}]

SIC Curve with statistical tests

SIC with centerline

Here we measure the centerline proposed by Harding, Leblanc, Goulet, & Cousineau (submitted). It is the median of all experimental conditions pooled together.

centerline=Median[Flatten[{{RThh,RThl,RTlh,RTll}}]];

Now we can measure and plot the SIC (along with the SIC centerline) with:

SIC=(shh+sll)-(shl+slh);

ListPlot[SIC,Joined->True,PlotRange->{All,{-1,1}},AxesLabel->{"t ","SIC(t)"},GridLines->{{{centerline,{Red,Dashed}}}}]

Kolmogorov-Smirnov test

Here we use an adapted Kolmogorov-Smirnov test (Houpt & Townsend, 2010) to see whether the SIC differs from 0. As we already have the SIC, which is shown above, all we need is the Harmonic Mean.

hMean=HarmonicMean[{Length[RThh],Length[RTlh],Length[RThl],Length[RTll]}];

Here we test the positive portion of the curve

(SuperPlus[D])=Max[SIC]
(SuperPlus[P])=Exp[-2(hMean/4)((((SuperPlus[D])))^2)]

Here we test the negative portion of the curve

(SuperMinus[D])=Min[SIC]
(SuperMinus[P])=Exp[-2(hMean/4)((((SuperMinus[D])))^2)]

Capacity

Here we test the capacity of our cognitive process. First, we must gather RT when a sub-process is working alone. In these cells, “RTh1”, “RTh2” are the individual RT for both sub-processes working alone in their H condition and “RThh” is the systems’s RT when both sub-processes are operating together.

RTh1=Table[Accumulate[RandomVariate[𝒟[(Subscript[μ, h1]),σ],{k}]]//Max,{n}];
RTh2=Table[Accumulate[RandomVariate[𝒟[(Subscript[μ, h2]),σ],{k}]]//Max,{n}];

time1=RandomVariate[𝒟[(Subscript[μ, h1]),σ],{n,k}];
time2=RandomVariate[𝒟[(Subscript[μ, h2]),σ],{n,k}];
RThh=MapThread[Coactive[#1,#2,k]&,{time1,time2}];

Finally, using the cumulative hazard function of each of the three working conditions (RTh1, RTh2, RThh) we are able to calculate and plot the capacity curve using CapOR as introduced in Townsend & Nozawa (1995) and CapOR1, a variant of CapOR that is shown in Houpt, Blaha, McIntire, Havig, & Townsend (2014). We use CapOR as we are working in a self-terminating environment.

hi=Max[Flatten[{{RThh,RTh1,RTh2}}]];
bins={0,hi,1};

Hab=HistogramList[RThh,bins,"CHF"][[2]];
Ha=HistogramList[RTh1,bins,"CHF"][[2]];
Hb=HistogramList[RTh2,bins,"CHF"][[2]];

Off[Power::infy];
Off[Infinity::indet];

CapOR=(Hab/((Ha+Hb)));
CapOR1=Hab-(Ha+Hb);


ListPlot[CapOR,Joined->True,PlotRange->{All,All},AxesLabel->{"t ","C_(OR)"}]
ListPlot[CapOR1,Joined->True,PlotRange->{All,All},AxesLabel->{"t ","C_(OR)"}]

The “Off” command above ensures that the calculations keep going regardless of any (1/0)calculations that may be present.