% caculate exchange reactions
exchange_rxns = [];
modelDir = 'data/gap-filling/iterative'
habitats = {'coral','seawater'};
methods = {'carveme','gapseq','kbase','comsensus'}

for i = 1:numel(habitats)
    for j = 1:numel(methods) 
        model_workspace = fullfile(modelDir, habitats{i}, [methods{j}, '.mat'])
        load(model_workspace,'GF')
            for k = 1:numel(GF)
                [selExc, selUpt] = findExcRxns(GF{k});
                if i == 1
                    exchange_rxns(k,j) = sum(selExc);
                else 
                   exchange_rxns(k,j+4) = sum(selExc);
                end
            end
    end
end

exchange_rxns = array2table(exchange_rxns);
exchange_rxns.Properties.VariableNames = {'c_carve','c_gapseq','c_kbase','c_merge','s_carve','s_gapseq','s_kbase','s_merge'};
writetable(exchange_rxns,'figures/exchange_rxns.csv');
