function example()
%% set a set of parameters
% generate the training data
nProts = 3;
nTimepoints = 10;
nConditions = 5;
data.data = getData(nProts,nTimepoints,nConditions);
% the measurement time points
data.t = linspace(0,1,nTimepoints)';
% the names of the all proteins involved in the signaling pathway
data.protNames = {'X1','X2','X3'};
% the names of the different stimuli perturbations (training data set)
data.conditionNames = {'WT','KO1','KO2','KO3','KO4'};

% visualise the results?
data.viz = false;

%% predict behaviour under unforeseen conditions
% just predict the network dynamics
data.infer = false;
% the connectivity matrix of the signaling pathway (prot_j -> prot_i := (i,j))      
% the column labels are: {'X1','X2','X3'}                  
data.networkMatrix = [0 1 1; % X1
                      1 1 0; % X2
                      1 1 1]; % X3
% generate the test data
data.testData = getData(nProts,nTimepoints,nConditions);
% names of the different stimuli perturbations (training data set)
data.conditionNamesTest = {'WT (test)','WT2 (test)','WT3 (test)','WT4 (test)','WT5 (test)'};
% estimate the trajectories
result = Sorad(data);
result.trajectories

%% calculate the marginal likelihoods of different models
% let us infer the network topology
data.infer = true;
% the maximum number of regulators (<= number of proteins)
data.maxRegulators = 3;
% calculate the marginal likelihoods
result = Sorad(data);
% the models, i.e., the regulatory proteins
result.subPowerset
% marginal likelihoods, in the rows are the models and in the
% columns are the target proteins
result.margLikelihoods

end

function data = getData(nProt,nTimepoints,nConditions)
% generate the data
data = cell(nProt,nConditions);
for protIdx = 1:nProt
    for condIdx = 1:nConditions
        data{protIdx,condIdx} = rand(nTimepoints,1);
    end
end
end