diff --git a/IQM Tools/IQMlite/Contents.m b/IQM Tools/IQMlite/Contents.m
index 4e9614b26812ad9ec926e882415458ed674a093f..ef8c5395750ba5d1519ae71d1dfe96de488bb59f 100644
--- a/IQM Tools/IQMlite/Contents.m
+++ b/IQM Tools/IQMlite/Contents.m
@@ -1,5 +1,5 @@
% IQM Tools Lite
-% Version 1.2.1 (R2015B) 30.04.2016
+% Version 1.2.2 (R2015B) 02.01.2017
%
% IQMlite
% =======
@@ -171,6 +171,7 @@
% IQMplot - IQMplot - plots given data.
% IQMplot2 - IQMplot2 - plots bar diagrams for given data.
% IQMplotCatCat - This function plots a "bubble" plot for categorical data to assess correlations
+% IQMplotKM - Kaplan-Meier plot
%
% functions
% These functions can be used on IQMmodels but on command line as well. When used in models then these are also working when converting the IQMmodel to C-code
@@ -301,7 +302,16 @@
% kin_uncomp_inihib_irr - Uncompetitive inhibition (irreversible) kinetics
% kin_uncomp_inihib_rev - Uncompetitive inhibition (reversible) kinetics
% kin_uni_uni_rev - Uni uni reversible kinetics
-%
+%
+% survival
+% IQMcoxExt - Parameter estimation for the extended Cox model
+% IQMcoxPH - Parameter estimation for the Cox proportional hazards model
+% IQMcoxStratPH - Parameter estimation for the Stratified Cox proportional hazards model
+% IQMhr - Hazard ratio estimate for Cox models
+% IQMkm - Kaplan-Meier curves estimation and plot (see also IQMplotKM)
+% IQMphTest - Schoenfeld residuals test for the proportional hazards assumption
+% IQMxext - Extension to time varying covariates (input to IQMcoxExt)
+%
% Auxiliary
% =========
% getDefaultIntegratorOptionsIQM - Inside this function/file default values for integrator options are set. This function also retrieves them
diff --git a/IQM Tools/IQMlite/auxiliary/compliance/IQMoutputTable.m b/IQM Tools/IQMlite/auxiliary/compliance/IQMoutputTable.m
new file mode 100644
index 0000000000000000000000000000000000000000..4fbd284b7b9c601306b385c84db777e09822001d
--- /dev/null
+++ b/IQM Tools/IQMlite/auxiliary/compliance/IQMoutputTable.m
@@ -0,0 +1,61 @@
+function [text] = IQMoutputTable(xtable,xtitle,xfooter,filename)
+% This function allows to convert a MATLAB table to a report table
+% formatted to be compatible with IQReport.
+%
+% [SYNTAX]
+% [text] = IQMoutputTable(xtable,xtitle,xfooter,filename)
+%
+% [INPUT]
+% xtable: MATLAB table
+% xtitle: String with the table caption (default: '')
+% xfooter: String with the table footer (default: '')
+% filename: If a filename is provided then the table is exported to
+% filename.txt. (default: '');
+%
+% [OUTPUT]
+% The text is written to a file.
+
+% <<<COPYRIGHTSTATEMENT - IQM TOOLS LITE>>>
+
+% Handle variable input arguments
+if nargin<4,
+ error('Please provide all input arguments.');
+end
+
+% Check table
+if ~istable(xtable),
+ error('First input argument is not a MATLAB table.');
+end
+
+% Get info
+Ncols = size(xtable,2)+1;
+headerNames = xtable.Properties.VariableNames;
+
+% Generate table-content as cell-array
+tableCells = table2cell(xtable);
+
+% Add <TR> info
+tableRows = cell(size(tableCells,1),1);
+tableRows(1:end) = {'<TR>'};
+tableCells = [tableRows tableCells];
+
+% Add <TH> info
+tableCells = [ [{'<TH>'} headerNames]; tableCells];
+
+% Add <TT> info
+tableTitle = cell(1,Ncols);
+tableTitle{1} = '<TT>';
+tableTitle{2} = xtitle;
+tableCells = [tableTitle; tableCells];
+
+% Add <TF> info
+if ~isempty(xfooter),
+ tableFooter = cell(1,Ncols);
+ tableFooter{1} = '<TF>';
+ tableFooter{2} = xfooter;
+ tableCells = [tableCells; tableFooter];
+end
+
+% Convert to report table
+IQMconvertCellTable2ReportTable(tableCells,'report',filename);
+
diff --git a/IQM Tools/IQMlite/auxiliary/graphics/ps/gs/fonts/EMPTY.txt b/IQM Tools/IQMlite/auxiliary/graphics/ps/gs/fonts/EMPTY.txt
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/IQM Tools/IQMlite/auxiliary/graphics/ps/gs/lib/EMPTY.txt b/IQM Tools/IQMlite/auxiliary/graphics/ps/gs/lib/EMPTY.txt
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/IQM Tools/IQMlite/auxiliary/output/IQMgetcolors.m b/IQM Tools/IQMlite/auxiliary/output/IQMgetcolors.m
index c2923b4c9b516823c98f5dad3dee61483e113f55..c9e5d2051d705c244c85e8be576229b3e9c423be 100644
--- a/IQM Tools/IQMlite/auxiliary/output/IQMgetcolors.m
+++ b/IQM Tools/IQMlite/auxiliary/output/IQMgetcolors.m
@@ -86,7 +86,7 @@ colors = [
% 0.85 0.7 1
% 1 0.6 0.78];
-lines = {'o-','x-','+-','*-','s-','d-','v-','^-','<-','>-','p-','h-', 'o--','x--','+--','*--','s--','d--','v--','^--','<--','>--','p--','h--', 'o-.','x-.','+-.','*-.','s-.','d-.','v-.','^-.','<-.','>-.','p-.','h-.', 'o:','x:','+:','*:','s:','d:','v:','^:','<:','>:','p:','h:' };
+lines = {'-' '--' ':' '-.' 'o-','x-','+-','*-','s-','d-','v-','^-','<-','>-','p-','h-', 'o--','x--','+--','*--','s--','d--','v--','^--','<--','>--','p--','h--', 'o-.','x-.','+-.','*-.','s-.','d-.','v-.','^-.','<-.','>-.','p-.','h-.', 'o:','x:','+:','*:','s:','d:','v:','^:','<:','>:','p:','h:' };
dots = {'o','x','+','*','s','d','v','^','<','>','p','h'};
diff --git a/IQM Tools/IQMlite/auxiliary/usernameIQM.m b/IQM Tools/IQMlite/auxiliary/usernameIQM.m
index fdc03d5e019c617e378aae9f97b376c2ed77a165..d1d87296b2180e6b80ae0cb81b9697dbbf7acde7 100644
--- a/IQM Tools/IQMlite/auxiliary/usernameIQM.m
+++ b/IQM Tools/IQMlite/auxiliary/usernameIQM.m
@@ -1,4 +1,4 @@
-function [username] = usernameIQM(text)
+function [username] = usernameIQM()
% usernameIQM: Get the name of the current user
% <<<COPYRIGHTSTATEMENT - IQM TOOLS LITE>>>
diff --git a/IQM Tools/IQMlite/classeshandling/modelhandling/auxiliary/text2IQMmodel/convertTextBCToModelIQM.m b/IQM Tools/IQMlite/classeshandling/modelhandling/auxiliary/text2IQMmodel/convertTextBCToModelIQM.m
index 54955e19de5775ffea240e712b56a3b4bff050c2..e7f0b4a5f912030a05a8fbd64c6a71c385463d06 100644
--- a/IQM Tools/IQMlite/classeshandling/modelhandling/auxiliary/text2IQMmodel/convertTextBCToModelIQM.m
+++ b/IQM Tools/IQMlite/classeshandling/modelhandling/auxiliary/text2IQMmodel/convertTextBCToModelIQM.m
@@ -698,7 +698,7 @@ for k=1:length(stateConstraintInfo),
highbound = stateConstraintInfo(k).highbound;
ODE = stateConstraintInfo(k).ODE;
% add a factor "constraints_factor_statename" to the ODE for the state
- IQMstructure.states(k).ODE = ['constraints_factor_' statename ' * (' ODE ')'];
+ IQMstructure.states(stateindex).ODE = ['constraints_factor_' statename ' * (' ODE ')'];
% define the "constraints_factor_statename" in the reactions
IQMstructure.reactions(end+1).name = ['constraints_factor_' statename];
% construct the piecewise expression
diff --git a/IQM Tools/IQMlite/classeshandling/modelhandling/auxiliary/text2IQMmodel/convertTextToModelIQM.m b/IQM Tools/IQMlite/classeshandling/modelhandling/auxiliary/text2IQMmodel/convertTextToModelIQM.m
index 74dbccf81cbea321789a2a7389543bcc3bbad148..aa5dbe8bc52e35232bdfa6f5c60eb51e31bd9976 100644
--- a/IQM Tools/IQMlite/classeshandling/modelhandling/auxiliary/text2IQMmodel/convertTextToModelIQM.m
+++ b/IQM Tools/IQMlite/classeshandling/modelhandling/auxiliary/text2IQMmodel/convertTextToModelIQM.m
@@ -132,7 +132,7 @@ for k=1:length(stateConstraintInfo),
highbound = stateConstraintInfo(k).highbound;
ODE = stateConstraintInfo(k).ODE;
% add a factor "constraints_factor_statename" to the ODE for the state
- IQMstructure.states(k).ODE = ['constraints_factor_' statename ' * (' ODE ')'];
+ IQMstructure.states(stateindex).ODE = ['constraints_factor_' statename ' * (' ODE ')'];
% define the "constraints_factor_statename" in the reactions
IQMstructure.reactions(end+1).name = ['constraints_factor_' statename];
% construct the piecewise expression
diff --git a/IQM Tools/IQMlite/release.txt b/IQM Tools/IQMlite/release.txt
index 2fbaf0c5ef2f71678c6d0cf65d8a7c5cac1ccaf4..4144d8a335062e61f86d1351a4a7e3cbbaa8d779 100644
--- a/IQM Tools/IQMlite/release.txt
+++ b/IQM Tools/IQMlite/release.txt
@@ -1,6 +1,12 @@
Release Notes IQM Tools Lite
============================
+Version 1.2.2 (02.01.2017)
+--------------------------
+- Minor big fixes
+- Improved graphical output
+- Survival analysis functions
+
Version 1.2.1 (30.04.2016)
------------------------
- Minor big fixes
diff --git a/IQM Tools/IQMlite/tools/dataset/IQMloadCSVdataset.m b/IQM Tools/IQMlite/tools/dataset/IQMloadCSVdataset.m
index 9bdc8430a41cdeab2bb4bb562fb54ec071a6d7ce..016ff07ee416ecaeb16aa3d7688f744481fcff92 100644
--- a/IQM Tools/IQMlite/tools/dataset/IQMloadCSVdataset.m
+++ b/IQM Tools/IQMlite/tools/dataset/IQMloadCSVdataset.m
@@ -190,11 +190,25 @@ if FLAG_HANDLE_SPECIAL_COLUMNS,
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
- % NaN values in DURATION exchange to 0
+ % If DATEDAY is of class datetime it should be converted to cellarray
+ % of strings
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
try
- data.DURATION(isnan(data.DURATION)) = 0;
+ if strcmp(class(data.DATEDAY),'datetime'),
+ data.DATEDAY = cellstr(data.DATEDAY);
+ end
+ end
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % If VISNAME is numeric, convert to cell with empty
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ try
+ if isnumeric(data.VISNAME),
+ data.VISNAME = cell(height(data),1);
+ data.VISNAME(1:end) = {''};
+ end
end
+
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
diff --git a/IQM Tools/IQMlite/tools/plots/IQMplotCatCat.m b/IQM Tools/IQMlite/tools/plots/IQMplotCatCat.m
index db88c468ef09a4cb3866f18482aa0fdb4a3515bc..7f35f3058d2da729bc2a0b99bf30d1d2db7601c8 100644
--- a/IQM Tools/IQMlite/tools/plots/IQMplotCatCat.m
+++ b/IQM Tools/IQMlite/tools/plots/IQMplotCatCat.m
@@ -1,10 +1,16 @@
-function [] = IQMplotCatCat(data,catNames)
+function [] = IQMplotCatCat(data,catNames,compareNames,options)
% This function plots a "bubble" plot for categorical data. This allows to
% at least graphically explore correlations between different categorical
% covariates. Areas of circles are proportional to the fraction.
%
+% There are 2 cases. Either all categorical covariates are compared to each
+% other (only catNames provided). Or categorical covariates are compared to
+% categorical covariates, defined by compareNames.
+%
% [SYNTAX]
% [] = IQMplotCatCat(data,catNames)
+% [] = IQMplotCatCat(data,catNames,compareNames)
+% [] = IQMplotCatCat(data,catNames,compareNames,options)
%
% [INPUT]
% data: Matlab dataset. Each column corresponds to a variable
@@ -12,7 +18,21 @@ function [] = IQMplotCatCat(data,catNames)
% defined in "catNames" need to be present in
% the dataset and contain numerical categorical data.
% catNames: Cell-array with names of categorical variables
-% options: MATLAB structure with optional arguments
+% compareNames: Cell-array with names of categorical variables to compare
+% with the categorical covariates in catNames
+% options: Matlab structure with additional options
+% options.percent: Show percentages instead of numbers. Percentages
+% are calculated relative to the number in the bins
+% of the catNames and only displayed if compared to
+% compareNames. =0 do show numbers and percent. =1 show
+% percent (default) =2 show numbers and percent.
+% options.bubbleSize: Factor to change size of bubbles (default: 500)
+% options.fontSizeText: Fontsize for the number and percent text
+% (default: 10)
+% options.hideFirstRow: Hides first row with bars in case
+% compareNames is used. =0: not hide (default),
+% =1: hide (default)
+% options.color: Color of bubbles default: 0.4*[1 1 1]
%
% [OUTPUT]
% Plot
@@ -24,74 +44,253 @@ if ~istable(data),
error('Input data needs to be provided as MATLAB table.');
end
-% Define plot color
-color = 0.4*[1 1 1];
-
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-% Check if dataset contains defined columns
+% Check which case it is and handle each case completely separately
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-for k=1:length(catNames),
- try
- data.(catNames{k});
- catch
- error(sprintf('Please check if "%s" is a column in the dataset!',catNames{k}));
- end
+if nargin==2,
+ compareNames = {};
+ options = [];
end
-
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-% Keep only selected columns
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-data = data(:,catNames);
-
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-% If is table then convert to double
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-if istable(data),
- data = table2array(data);
+if nargin==3,
+ options = [];
end
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-% Keep only selected columns
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-
-% % Default settings
-% % Flag to take log of values
-% LogFlag = 0;
-% % Coloring
-% CorrThres = 0.3; % Pearson correlation coefficient threshold for different background color
-% AxisColor = [1 .2 .2]; %color of axis when
-%
-% % Get optional values if defined
-% try names = OPTIONS.names; catch, end;
-% try LogFlag = OPTIONS.LogFlag; catch, end;
-% try CorrThres = OPTIONS.CorrThres; catch, end;
-% try AxisColor = OPTIONS.AxisColor; catch, end;
+try percentFlag = options.percent; catch, percentFlag = 1; end
+try bubbleSize = options.bubbleSize; catch, bubbleSize = 500; end
+try fontSizeText = options.fontSizeText; catch, fontSizeText = 10; end
+try hideFirstRow = options.hideFirstRow; catch, hideFirstRow = 0; end
+try color = options.color; catch, color = 0.4*[1 1 1]; end
% Subindex properties
-Spacing = 0;
-Padding = 0;
+Spacing = 0.0005;
+Padding = 0.001;
Margin = .1;
% Define fontsize
if length(catNames)<=6,
FONTSIZE = 10;
else
- FONTSIZE = 8;
+ FONTSIZE = 10;
end
-clf;
-n = length(catNames);
-for ir=1:n
- ystr = catNames{ir};
- y = data(:,ir);
- for ic=1:ir
- xstr = catNames{ic};
- x = data(:,ic);
- ip = (ir-1)*n+ic;
- subaxis(n,n,ip,'Spacing',Spacing,'Padding',Padding,'Margin',Margin);
- if ir==ic % plot "histogram"
-
+if isempty(compareNames),
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % Handle first case (compare all catNames with each other)
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % Check if dataset contains defined columns
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ for k=1:length(catNames),
+ try
+ data.(catNames{k});
+ catch
+ error(sprintf('Please check if "%s" is a column in the dataset!',catNames{k}));
+ end
+ end
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % Keep only selected columns
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ data = data(:,catNames);
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % If is table then convert to double
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ if istable(data),
+ data = table2array(data);
+ end
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % Keep only selected columns
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+ clf;
+ n = length(catNames);
+ for ir=1:n
+ ystr = catNames{ir};
+ y = data(:,ir);
+ for ic=1:ir
+ xstr = catNames{ic};
+ x = data(:,ic);
+ ip = (ir-1)*n+ic;
+ subaxis(n,n,ip,'Spacing',Spacing,'Padding',Padding,'Margin',Margin);
+ if ir==ic % plot "histogram"
+
+ xu = unique(x);
+ xu = xu(~isnan(xu));
+ b = zeros(size(xu));
+ nu = length(xu);
+ for i=1:nu
+ b(i) = sum(x==xu(i));
+ end
+ h = bar(1:nu,b,0.5);
+
+ % Write out numbers
+ for knr=1:length(b),
+ aaa = get(gca,'YLim');
+ text(knr,b(knr)+max(aaa)*0.02,sprintf('%d',b(knr)),'VerticalAlign','bottom','HorizontalAlign','center','FontSize',fontSizeText+2);
+ end
+
+ set(h,'FaceColor',color)
+ set(gca,'XLim',[0.5 nu+.5]);
+ title(xstr,'Interpreter','none')
+ set(gca,'YTick',[]);
+ set(gca,'XTick',[]);
+ set(gca,'YLim',[0 max(b)*1.3]);
+ set(gca,'FontSize',fontSizeText);
+
+ if ir==1,
+ ylabel('#','FontSize',8);
+ end
+
+ if ir==n
+ xlabel(xstr,'Interpreter','none','FontSize',fontSizeText);
+ set(gca,'XTick',1:length(xu))
+ set(gca,'XTickLabel',xu);
+ end
+
+ else % plot bubble plot
+ % Exchange x categories for numbers 1:N
+ xu = unique(x);
+ iu = {};
+ for k=1:length(xu),
+ iu{k} = find(x==xu(k));
+ end
+ xn = NaN(size(x));
+ for k=1:length(iu),
+ xn(iu{k}) = k;
+ end
+
+ % Exchange y categories for numbers 1:N
+ yu = unique(y);
+ iu = {};
+ for k=1:length(yu),
+ iu{k} = find(y==yu(k));
+ end
+ yn = NaN(size(y));
+ for k=1:length(iu),
+ yn(iu{k}) = k;
+ end
+
+ % Find unique pairs of xn and yn (these will be plotted in the
+ % bubble plot
+ [pairs,~,b] = unique([xn yn],'rows');
+
+ % Determine number of rows in [xn yn] with these pairs
+ nr_pairs = NaN(size(pairs,1),1);
+
+ for k=1:size(pairs,1),
+ nr_pairs(k) = sum(b==k);
+ end
+
+ % Determine relative number of pairs in fraction
+ nr_pairs_fraction = nr_pairs/sum(nr_pairs);
+
+ % Determine size based on fraction
+ use_size = bubbleSize*nr_pairs_fraction;
+
+ % Determine default color
+ color_use = color(ones(1,length(nr_pairs)),:);
+
+ % Plot
+ scatter(pairs(:,1),pairs(:,2), use_size, color_use, 'filled', 'MarkerEdgeColor', 'none');
+
+ % Annotate
+ set(gca,'XLim',[min(xn)-0.5 max(xn)+0.5])
+ set(gca,'YLim',[min(yn)-0.5 max(yn)+0.5])
+
+ % Write out percentage as text
+ for knr=1:size(pairs,1),
+ text(pairs(knr,1)+0.025,pairs(knr,2)+0.2,sprintf('%d',nr_pairs(knr)),'FontSize',FONTSIZE);
+ end
+
+ set(gca,'XTick',1:length(xu))
+ set(gca,'YTick',1:length(yu))
+ set(gca,'FontSize',fontSizeText);
+
+ if ic==1
+ ylabel(ystr,'Interpreter','none','FontSize',fontSizeText);
+ set(gca,'YTickLabel',yu);
+ else
+ set(gca,'YTickLabel',[]);
+ end
+
+ if ir==n
+ xlabel(xstr,'Interpreter','none','FontSize',fontSizeText);
+ set(gca,'XTickLabel',xu);
+ else
+ set(gca,'XTickLabel',[]);
+ end
+
+ grid on;
+
+ end
+ end
+ end
+
+else
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % Handle second case (compare all catNames with compareNames)
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % Check if dataset contains defined columns
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ for k=1:length(catNames),
+ try
+ data.(catNames{k});
+ catch
+ error(sprintf('Please check if "%s" is a column in the dataset!',catNames{k}));
+ end
+ end
+
+ for k=1:length(compareNames),
+ try
+ data.(compareNames{k});
+ catch
+ error(sprintf('Please check if "%s" is a column in the dataset!',compareNames{k}));
+ end
+ end
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % Keep only selected columns
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ dataCompare = data(:,compareNames);
+ data = data(:,catNames);
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % If is table then convert to double
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ if istable(data),
+ data = table2array(data);
+ end
+ if istable(dataCompare),
+ dataCompare = table2array(dataCompare);
+ end
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % Options / fontsizes
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+ clf;
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % First row is frequencies of catNames
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ if hideFirstRow==0,
+ ncols = length(catNames);
+ nrows = 1+length(compareNames);
+ for ir=1:ncols
+ ystr = catNames{ir};
+ y = data(:,ir);
+ xstr = catNames{ir};
+ x = data(:,ir);
+ ip = ir;
+ subaxis(nrows,ncols,ip,'Spacing',Spacing,'Padding',Padding,'Margin',Margin);
+ % Plot frequencies
xu = unique(x);
xu = xu(~isnan(xu));
b = zeros(size(xu));
@@ -100,24 +299,54 @@ for ir=1:n
b(i) = sum(x==xu(i));
end
h = bar(1:nu,b,0.5);
- set(h,'FaceColor',0.4*[1 1 1])
+ grid on
+ set(gca,'FontSize',fontSizeText);
+
+ % Write out numbers
+ for knr=1:length(b),
+ aaa = get(gca,'YLim');
+ text(knr,b(knr)+max(aaa)*0.02,sprintf('%d',b(knr)),'VerticalAlign','bottom','HorizontalAlign','center','FontSize',fontSizeText+2);
+ end
+
+ set(h,'FaceColor',color)
set(gca,'XLim',[0.5 nu+.5]);
title(xstr,'Interpreter','none')
- set(gca,'YTick',[]);
+% set(gca,'YTick',[]);
+ if ir>1,
+ set(gca,'YTickLabel','');
+ end
set(gca,'XTick',[]);
- set(gca,'YLim',[0 length(x)]);
+ set(gca,'YLim',[0 max(b)*1.3]);
if ir==1,
ylabel('#','FontSize',8);
end
+ end
+ offset = 1;
+ else
+ ncols = length(catNames);
+ nrows = length(compareNames);
+ offset = 0;
+ end
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % Next rows is correlation of catNames with compareNames
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % Cycle through compareNames
+ for kcompare=1:length(compareNames),
+ ir = kcompare+offset;
+ ystr = compareNames{kcompare};
+ y = dataCompare(:,kcompare);
+
+ % Cycle through covNames
+ for kcov=1:length(catNames),
+ ic = kcov;
+ xstr = catNames{ic};
+ x = data(:,ic);
+ ip = (ir-1)*ncols+ic;
+ subaxis(nrows,ncols,ip,'Spacing',Spacing,'Padding',Padding,'Margin',Margin);
- if ir==n
- xlabel(xstr,'Interpreter','none','FontSize',FONTSIZE);
- set(gca,'XTick',1:length(xu))
- set(gca,'XTickLabel',xu);
- end
-
- else % plot bubble plot
+ % plot bubble plot
% Exchange x categories for numbers 1:N
xu = unique(x);
iu = {};
@@ -128,7 +357,7 @@ for ir=1:n
for k=1:length(iu),
xn(iu{k}) = k;
end
-
+
% Exchange y categories for numbers 1:N
yu = unique(y);
iu = {};
@@ -138,7 +367,7 @@ for ir=1:n
yn = NaN(size(y));
for k=1:length(iu),
yn(iu{k}) = k;
- end
+ end
% Find unique pairs of xn and yn (these will be plotted in the
% bubble plot
@@ -155,31 +384,38 @@ for ir=1:n
nr_pairs_fraction = nr_pairs/sum(nr_pairs);
% Determine size based on fraction
- use_size = 500*nr_pairs_fraction;
+ use_size = bubbleSize*nr_pairs_fraction;
% Determine default color
color_use = color(ones(1,length(nr_pairs)),:);
% Plot
scatter(pairs(:,1),pairs(:,2), use_size, color_use, 'filled', 'MarkerEdgeColor', 'none');
-
+
% Annotate
set(gca,'XLim',[min(xn)-0.5 max(xn)+0.5])
set(gca,'YLim',[min(yn)-0.5 max(yn)+0.5])
+ % Write out numbers as text
+ for knr=1:size(pairs,1),
+ textShow = sprintf('%d',nr_pairs(knr));
+ text(pairs(knr,1)+0.025,pairs(knr,2)+0.2,textShow,'FontSize',FONTSIZE);
+ end
+
set(gca,'XTick',1:length(xu))
- set(gca,'YTick',1:length(yu))
-
+ set(gca,'YTick',1:length(yu))
+ set(gca,'FontSize',fontSizeText);
+
if ic==1
- ylabel(ystr,'Interpreter','none','FontSize',FONTSIZE);
+ ylabel(ystr,'Interpreter','none','FontSize',10);
set(gca,'YTickLabel',yu);
else
set(gca,'YTickLabel',[]);
- end
-
- if ir==n
- xlabel(xstr,'Interpreter','none','FontSize',FONTSIZE);
- set(gca,'XTickLabel',xu);
+ end
+
+ if ir==nrows
+ xlabel(xstr,'Interpreter','none','FontSize',fontSizeText);
+ set(gca,'XTickLabel',xu);
else
set(gca,'XTickLabel',[]);
end
@@ -188,4 +424,21 @@ for ir=1:n
end
end
+
+ % Adjust YLim on each row to same values
+ if length(catNames) > 1,
+ for krow=1:length(compareNames)+1,
+ YLimALL = [];
+ for kcol=1:length(catNames),
+ subaxis(length(compareNames)+1,length(catNames),(krow-1)*length(catNames)+kcol);
+ YLimALL = [YLimALL; get(gca,'YLim')];
+ end
+ YLimMin = min(YLimALL);
+ YLimMax = max(YLimALL);
+ for kcol=1:length(catNames),
+ subaxis(length(compareNames)+1,length(catNames),(krow-1)*length(catNames)+kcol);
+ set(gca,'YLim',[YLimMin(1) YLimMax(2)]);
+ end
+ end
+ end
end
\ No newline at end of file
diff --git a/IQM Tools/IQMlite/tools/plots/IQMplotCovarianceCat.m b/IQM Tools/IQMlite/tools/plots/IQMplotCovarianceCat.m
index 9d6207b393a5d5e3e91a1ee22203379ee941f996..b654e36c88015bba572f5ddf4d9b181325b3dcd0 100644
--- a/IQM Tools/IQMlite/tools/plots/IQMplotCovarianceCat.m
+++ b/IQM Tools/IQMlite/tools/plots/IQMplotCovarianceCat.m
@@ -1,7 +1,7 @@
function [] = IQMplotCovarianceCat(data,contNames,catNames,options)
-% This function plots the covariance relationship between a list of
+% This function plots the covariance relationship between a list of
% continuous variables (contNames) and a list of categorical variables
-% (catNames), passed in "data".
+% (catNames), passed in "data".
%
% [SYNTAX]
% [] = IQMplotCovarianceCat(data,contNames,catNames)
@@ -11,13 +11,15 @@ function [] = IQMplotCovarianceCat(data,contNames,catNames,options)
% data: Matlab dataset. Each column corresponds to a variable
% and each row to a sample. The columns with the names
% defined in "contNames" and "catNames" need to be present in
-% the dataset.
+% the dataset.
% contNames: Cell-array with names of continuous variables
% catNames: Cell-array with names of categorical variables
% options: MATLAB structure with optional arguments
%
% options.LogFlag: =1 do log transform the variables,
% =0 do not transform (default: 0)
+% options.fontSizeText: Fontsize for the number and percent text
+% (default: 10)
%
% [OUTPUT]
% Plot
@@ -26,21 +28,22 @@ function [] = IQMplotCovarianceCat(data,contNames,catNames,options)
LogFlag = 0;
try LogFlag = options.LogFlag; catch, end
+try fontSizeText = options.fontSizeText; catch, fontSizeText = 10; end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Check if dataset contains defined columns
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for k=1:length(contNames),
- try
- data.(contNames{k});
- catch,
+ try
+ data.(contNames{k});
+ catch,
error(sprintf('Please check if "%s" is a column in the dataset!',contNames{k}));
end
end
for k=1:length(catNames),
- try
- data.(catNames{k});
- catch,
+ try
+ data.(catNames{k});
+ catch,
error(sprintf('Please check if "%s" is a column in the dataset!',catNames{k}));
end
end
@@ -57,7 +60,7 @@ Margin = .1;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
ncat = length(catNames); %rows are the categorical covariates
ncts = length(contNames);
-for icat=1:ncat
+for icat=1:ncat
xstr = catNames{icat};
x = data.(xstr);
@@ -70,6 +73,7 @@ for icat=1:ncat
ip=icat;
subaxis(ncts+1,ncat,ip,'Spacing',Spacing,'Padding',Padding,'Margin',Margin);
+ grid on
xu = unique(x);
xu = xu(~isnan(xu));
@@ -79,17 +83,28 @@ for icat=1:ncat
b(i) = sum(x==xu(i));
end
h = bar(1:nu,b,0.5);
- set(h,'FaceColor',0.4*[1 1 1])
+ set(h,'FaceColor',0.6*[1 1 1])
set(gca,'XLim',[0.5 nu+.5]);
+ set(gca,'FontSize',fontSizeText);
+
+ grid on
+
+ % Write out numbers
+ for knr=1:length(b),
+ aaa = get(gca,'YLim');
+ text(knr,b(knr)+max(aaa)*0.02,sprintf('%d',b(knr)),'VerticalAlign','bottom','HorizontalAlign','center','FontSize',fontSizeText+2);
+ end
+
title(xstr,'Interpreter','none')
if icat==1
ylabel('#','Interpreter','none');
else
- set(gca,'YTick',[]);
+ % set(gca,'YTick',[]);
+ set(gca,'YTickLabel','');
end
set(gca,'XTick',[]);
- set(gca,'YLim',[0 length(x)]);
+ set(gca,'YLim',[0 max(b)*1.3]);
for icts=1:ncts
ystr = contNames{icts};
@@ -112,7 +127,7 @@ for icat=1:ncat
M(x==xx(i),i) = y(x==xx(i));
end
OPTIONSbox.NumFlag = 0;
- OPTIONSbox.BoxColor = 0.4*[1 1 1];
+ OPTIONSbox.BoxColor = 0.6*[1 1 1];
OPTIONSbox.BoxWidth = .5;
OPTIONSbox.MedianWidth = .7;
OPTIONSbox.OutlierColor = 0.4*[1 1 1];
@@ -120,6 +135,8 @@ for icat=1:ncat
OPTIONSbox.MedianColor = [0 0 0];
plotboxIQM(M,1:length(xx),OPTIONSbox);
set(gca,'XLim',[.5 length(xx)+.5]);
+ grid on
+ set(gca,'FontSize',fontSizeText);
set(gca,'XTick',1:length(xx));
if icts<ncts
@@ -132,8 +149,11 @@ for icat=1:ncat
end
if icat==1
ylabel(ystr,'Interpreter','none');
+ else
+ % set(gca,'YTick',[]);
+ set(gca,'YTickLabel','');
end
- set(gca,'YTick',[]);
+ % set(gca,'YTick',[]);
if min(y)~=max(y),
set(gca,'YLim',[min(y) max(y)]);
else
@@ -141,6 +161,21 @@ for icat=1:ncat
end
end
end
- set(gca,'YTick',[]);
+ % set(gca,'YTick',[]);
+ end
+end
+
+% Adjust YLim on each row to same values
+for krow=1:ncts+1,
+ YLimALL = [];
+ for kcol=1:ncat,
+ subaxis(ncts+1,ncat,(krow-1)*ncat+kcol);
+ YLimALL = [YLimALL; get(gca,'YLim')];
end
+ YLimMin = min(YLimALL);
+ YLimMax = max(YLimALL);
+ for kcol=1:ncat,
+ subaxis(ncts+1,ncat,(krow-1)*ncat+kcol);
+ set(gca,'YLim',[min(YLimMin) max(YLimMax)]);
+ end
end
diff --git a/IQM Tools/IQMlite/tools/plots/IQMplotKM.m b/IQM Tools/IQMlite/tools/plots/IQMplotKM.m
new file mode 100644
index 0000000000000000000000000000000000000000..200b517a248d53976ccfad91737f7847eecf6709
--- /dev/null
+++ b/IQM Tools/IQMlite/tools/plots/IQMplotKM.m
@@ -0,0 +1,101 @@
+function [ST50_output_cellTable,ST50_output] = IQMplotKM(TIME,CENS,color,line,type,marksFlag,CIflag)
+% Plots a Kaplan Meier curve.
+%
+% [SYNTAX]
+% [ST50_output] = IQMplotKM(TIME,CENS)
+% [ST50_output] = IQMplotKM(TIME,CENS,color)
+% [ST50_output] = IQMplotKM(TIME,CENS,color,line)
+% [ST50_output] = IQMplotKM(TIME,CENS,color,line,type)
+% [ST50_output] = IQMplotKM(TIME,CENS,color,line,type,marksFlag)
+% [ST50_output] = IQMplotKM(TIME,CENS,color,line,type,marksFlag,CIflag)
+%
+% [INPUT]
+% TIME: vector with time information for time to event
+% CENS: vector with censoring information (0: uncensored, 1: right censored)
+% color: [r g b] color default: black
+% line: Line style (default: "-")
+% type: "Survivor" (default) or "cumulative hazard" or "cdf"
+% marksFlag: =1: censoring marks are plotted (default). =0 not plotted
+% CIflag: =1: plots 95% CI. =0 not plotted (default)
+%
+% [OUTPUT]
+% ST50_output: [50% survival time, 95% CI - lower and upper bound]
+
+% <<<COPYRIGHTSTATEMENT - IQM TOOLS LITE>>>
+
+if nargin<3,
+ color = [0 0 0];
+end
+if nargin<4,
+ line = '-';
+end
+if nargin<5,
+ type = 'survivor';
+end
+if nargin<6,
+ marksFlag = 1;
+end
+if nargin<7,
+ CIflag = 0;
+end
+
+% Remove NaN values from TIME and CENS (paired)
+ixNAN = find(isnan(TIME));
+TIME(ixNAN) = [];
+CENS(ixNAN) = [];
+
+% Plot KM curve
+[f,x,flo,fup] = ecdf(TIME,'censoring',CENS,'function',type);
+stairs(x,f,line,'LineWidth',1,'Color',color,'LineWidth',2); hold on;
+
+% Plot marks for censored data
+if marksFlag,
+ xi = TIME(logical(CENS));
+ [~,~,step] = histcounts(xi,x);
+ step(xi>max(TIME)) = length(x);
+ ixstep0 = find(step==0);
+ xi(ixstep0) = [];
+ step(ixstep0) = [];
+ fi = f(step);
+ plot(xi,fi, 's','MarkerEdgeColor','white','MarkerFaceColor','white','MarkerSize',8)
+ plot(xi,fi, 'b+','color',color,'MarkerSize',8,'LineWidth',1)
+end
+
+% Plot 95% confidence interval
+if CIflag,
+ IQMplotfill(x,flo,fup,color,0.1); hold on
+end
+
+% Annotate
+grid on
+set(gca,'FontSize',12);
+xlabel('Time','FontSize',14);
+ylabel('Probability','FontSize',14);
+set(gca,'YLim',[-0.05 1.05]);
+
+% Determine 50% survival time and 95% CI
+ix = find(f-0.5<0);
+if isempty(ix),
+ ST50 = NaN;
+else
+ ST50 = x(ix(1));
+end
+
+ix = find(flo-0.5<0);
+if isempty(ix),
+ ST50_CI_lo = NaN;
+else
+ ST50_CI_lo = x(ix(1));
+end
+
+ix = find(fup-0.5<0);
+if isempty(ix),
+ ST50_CI_up = NaN;
+else
+ ST50_CI_up = x(ix(1));
+end
+
+ST50_output = [ST50, ST50_CI_lo, ST50_CI_up];
+
+ST50_output_cellTable(1,1:3) = {'<TH>' '50% Survival Time' '95% CI'};
+ST50_output_cellTable(2,1:3) = {'<TR>' sprintf('%1.3g',ST50) sprintf('[%1.3g, %1.3g]',ST50_CI_lo,ST50_CI_up)};
diff --git a/IQM Tools/IQMlite/tools/plots/IQMplotPercentStratified.m b/IQM Tools/IQMlite/tools/plots/IQMplotPercentStratified.m
new file mode 100644
index 0000000000000000000000000000000000000000..4e51adcc17ad93076c9e98b7fc2c2f3b28d4cd3b
--- /dev/null
+++ b/IQM Tools/IQMlite/tools/plots/IQMplotPercentStratified.m
@@ -0,0 +1,200 @@
+function [yn_all_percent] = IQMplotPercentStratified(data,catNames,compareNames,options)
+% This function plots column plots / stacked bar plots for categorical outcomes (percent
+% of occurrence) stratified by categorical compare variables in
+% compareNames.
+%
+% Works for one catName and one compareName for now.
+%
+% [SYNTAX]
+% [] = IQMplotPercentStratified(data,catNames,compareNames)
+%
+% [INPUT]
+% data: Matlab dataset. Each column corresponds to a variable
+% and each row to a sample. The columns with the names
+% defined in "catNames" need to be present in
+% the dataset and contain numerical categorical data.
+% catNames: Cell-array with names of categorical variables
+% compareNames: Cell-array with names of categorical variables to compare
+% with the categorical covariates in catNames
+% options: Matlab structure with additional options
+% options.percent: Show percentages instead of numbers. Percentages
+% are calculated relative to the number in the bins
+% of the catNames and only displayed if compared to
+% compareNames. =0 do show numbers and percent. =1 show
+% percent (default) =2 show numbers and percent.
+% options.fontSizeText: Fontsize for the number and percent text
+% (default: 10)
+% options.color: Color of bubbles default: 0.4*[1 1 1]
+%
+% [OUTPUT]
+% Plot
+
+% <<<COPYRIGHTSTATEMENT - IQM TOOLS LITE>>>
+
+% Check if data is table
+if ~istable(data),
+ error('Input data needs to be provided as MATLAB table.');
+end
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Check which case it is and handle each case completely separately
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+if nargin==3,
+ options = [];
+end
+
+try fontSizeText = options.fontSizeText; catch, fontSizeText = 10; end
+try hideFirstRow = options.hideFirstRow; catch, hideFirstRow = 0; end
+try color = options.color; catch, color = 0.4*[1 1 1]; end
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Check if dataset contains defined columns
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+for k=1:length(catNames),
+ try
+ data.(catNames{k});
+ catch
+ error(sprintf('Please check if "%s" is a column in the dataset!',catNames{k}));
+ end
+end
+
+for k=1:length(compareNames),
+ try
+ data.(compareNames{k});
+ catch
+ error(sprintf('Please check if "%s" is a column in the dataset!',compareNames{k}));
+ end
+end
+
+if length(catNames)>1 || length(compareNames)>1,
+ error('Only single entry in catNames and compareNames allowed.');
+end
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Generate the plot
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Keep only selected columns
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+dataCompare = data(:,compareNames);
+data = data(:,catNames);
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% If is table then convert to double
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+if istable(data),
+ data = table2array(data);
+end
+if istable(dataCompare),
+ dataCompare = table2array(dataCompare);
+end
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Options / fontsizes
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+% Subindex properties
+Spacing = 0.00;
+Padding = 0.0005;
+Margin = .1;
+
+% Define fontsize
+if length(catNames)<=6,
+ FONTSIZE = 10;
+else
+ FONTSIZE = 8;
+end
+
+% Define ncols and
+ncols = length(catNames);
+nrows = length(compareNames);
+
+% Cycle through compareNames
+for kcompare=1:length(compareNames),
+ ir = kcompare;
+ ystr = compareNames{kcompare};
+ y = dataCompare(:,kcompare);
+
+ % Cycle through covNames
+ for kcov=1:length(catNames),
+ ic = kcov;
+ xstr = catNames{ic};
+ x = data(:,ic);
+ ip = (ir-1)*ncols+ic;
+% subaxis(nrows,ncols,ip,'Spacing',Spacing,'Padding',Padding,'Margin',Margin);
+
+ % Exchange x categories for numbers 1:N
+ xu = unique(x);
+ % Find positions where the unique items are listed (to obtain
+ % number 1...N) ... stored in xn
+ iu = {};
+ for k=1:length(xu),
+ iu{k} = find(x==xu(k));
+ end
+ xn = NaN(size(x));
+ for k=1:length(iu),
+ xn(iu{k}) = k;
+ end
+
+ % Exchange y categories for numbers 1:N
+ yu = unique(y);
+ iu = {};
+ for k=1:length(yu),
+ iu{k} = find(y==yu(k));
+ end
+ yn = NaN(size(y));
+ for k=1:length(iu),
+ yn(iu{k}) = k;
+ end
+
+ % Find numbers of unique yn in unique xn
+ yn_all_percent = [];
+ for k=1:max(xn),
+ % Get number of xn values in bin k
+ NR_xnk = length(xn(xn==k));
+ % Get number of different yn values in bin k
+ yn_all = [];
+ for k2=1:max(yn),
+ ynk = sum(xn==k & yn==k2);
+ yn_all(k2) = ynk;
+ end
+ % Determine percent
+ yn_all_percent = [yn_all_percent; 100*yn_all/NR_xnk];
+ end
+
+ % Plot
+ hp = bar(1:max(xn),yn_all_percent,'stacked');
+
+% hatch_combinations = {
+% 'fill' 0
+% 'single' 0
+% 'single' 45
+% 'single' 90
+% 'cross' 0
+% 'cross' 45
+% };
+%
+% for k=1:min(length(hp),6),
+% hatchfill2(hp(k),hatch_combinations{k,1},'HatchAngle',hatch_combinations{k,2});
+% end
+
+ % Annotate
+ if ic==1
+ ylabel('Percentage of x-bin','Interpreter','none','FontSize',12);
+ else
+ set(gca,'YTickLabel',[]);
+ end
+
+ if ir==nrows
+ xlabel(xstr,'Interpreter','none','FontSize',FONTSIZE);
+ set(gca,'XTickLabel',xu);
+ else
+ set(gca,'XTickLabel',[]);
+ end
+
+ set(gca,'YLim',[0 100]);
+ grid on;
+
+ end
+end
diff --git a/IQM Tools/IQMlite/tools/plots/IQMplotfacetgrid.m b/IQM Tools/IQMlite/tools/plots/IQMplotfacetgrid.m
index 4604d19843cea18e9db351e605a45d0a0ab1443f..0c5a8837b1b0ce705a2df0777559093d82075143 100644
--- a/IQM Tools/IQMlite/tools/plots/IQMplotfacetgrid.m
+++ b/IQM Tools/IQMlite/tools/plots/IQMplotfacetgrid.m
@@ -230,7 +230,7 @@ for kX=1:length(allGroupX),
if dataXYC.colorgroup(1) == -1,
datatemp.color = -1*ones(height(datatemp),1);
else
- datatemp.color = find(allGROUPcolor==dataXYC.colorgroup(1))*ones(length(datatemp),1);
+ datatemp.color = find(allGROUPcolor==dataXYC.colorgroup(1))*ones(height(datatemp),1);
end
else
datatemp.color = find(strcmp(allGROUPcolor,dataXYC.colorgroup{1}))*ones(height(datatemp),1);
diff --git a/IQM Tools/IQMlite/tools/plots/IQMplotpairwiseCorr.m b/IQM Tools/IQMlite/tools/plots/IQMplotpairwiseCorr.m
index 17eb0198510763e47dafe87861f5dbcf34da367f..ac2f85b795f90949c21a00a1ac4ec60665d18f3b 100644
--- a/IQM Tools/IQMlite/tools/plots/IQMplotpairwiseCorr.m
+++ b/IQM Tools/IQMlite/tools/plots/IQMplotpairwiseCorr.m
@@ -49,7 +49,7 @@ end
LogFlag = 0;
% Coloring
CorrThres = 0.3; % Pearson correlation coefficient threshold for different background color
-AxisColor = [1 .2 .2]; %color of axis when
+AxisColor = [1 .4 .4]; %color of axis when
% Get optional values if defined
try names = OPTIONS.names; catch, end;
@@ -111,7 +111,7 @@ for ir=1:n
if abs(rho)>CorrThres
optcorr.Color = AxisColor;
else
- optcorr.Color = 0.4*[1 1 1];
+ optcorr.Color = 0.6*[1 1 1];
end
optcorr.LineColor = [0 0 0];
optcorr.TitleType = 'none';
@@ -131,7 +131,7 @@ for ir=1:n
else
str = '|corr|<0.01';
end
- text(xt,yt,str,'Color',[0 0 0],'Hor','Center','Ver','Middle','FontWeight','Bold','Interpreter','none');
+ text(xt,yt,str,'Color',[0 0 0],'Hor','Right','Ver','Middle','FontWeight','Bold','Interpreter','none');
if min(xuse)~=max(xuse),
set(gca,'XLim',[min(xuse) max(xuse)]);
diff --git a/IQM Tools/IQMlite/tools/plots/auxiliary/plotgeneralIQM.m b/IQM Tools/IQMlite/tools/plots/auxiliary/plotgeneralIQM.m
index 191c97d48ba527e1c6d2d7e35dee8084256ac837..1d65309da87393609674fee1f15bee2d84fdedbe 100644
--- a/IQM Tools/IQMlite/tools/plots/auxiliary/plotgeneralIQM.m
+++ b/IQM Tools/IQMlite/tools/plots/auxiliary/plotgeneralIQM.m
@@ -363,13 +363,13 @@ if showloessline,
% Generate loess smoothed line
if logX==0 && logY==0,
- Ysmooth = smoothIQM(Xdata,Ydata,0.1,'loess');
+ Ysmooth = smoothIQM(Xdata,Ydata,0.5,'loess');
elseif logX==1 && logY==0,
- Ysmooth = smoothIQM(log(Xdata),Ydata,0.1,'loess');
+ Ysmooth = smoothIQM(log(Xdata),Ydata,0.5,'loess');
elseif logX==0 && logY==1,
- Ysmooth = exp(smoothIQM(Xdata,log(Ydata),0.1,'loess'));
+ Ysmooth = exp(smoothIQM(Xdata,log(Ydata),0.5,'loess'));
else
- Ysmooth = exp(smoothIQM(log(Xdata),log(Ydata),0.1,'loess'));
+ Ysmooth = exp(smoothIQM(log(Xdata),log(Ydata),0.5,'loess'));
end
% Handle color, etc.
if datakkk.color(1) == -1,
diff --git a/IQM Tools/IQMlite/tools/survival/IQMcoxExt.m b/IQM Tools/IQMlite/tools/survival/IQMcoxExt.m
new file mode 100644
index 0000000000000000000000000000000000000000..b7b77c491cd9302c405362f0b813cfa2e2fa53a0
--- /dev/null
+++ b/IQM Tools/IQMlite/tools/survival/IQMcoxExt.m
@@ -0,0 +1,202 @@
+function [stats,st] = IQMcoxExt(sID_ext,X_ext,d_ext,t_ext,e)
+% [stats,st] = IQMcoxExt(sID_ext,X_ext,d_ext,t_ext,e)
+%
+% Parameter estimation for the extended Cox model. This function uses as
+% input the output of IQMxext but with X_ext having new time-dependent
+% columns added by the user (eg. the product of a column of X_ext with
+% t_ext. Survival time ties are handled using Efron's method.
+%
+% Input
+% sID_ext: Extended subjects' IDs (nx1).
+% X_ext: Extended design matrix (nxp).
+% d_ext: Extended censorship status vector (nx1).
+% t_ext: Extended survival time vector (nx1).
+% e: Convergence epsilon (gradient's norm). Default 10^-3;
+%
+% Output
+% stats.Bhat: Estimated vector of the population regression parameters.
+% stats.CovBhat: Estimated covariance matrix of the population regression
+% parameters.
+% stats.llh: Values of the maximum log-likelihood across the optimization
+% process.
+% st: Termination state (1 for convergence and 0 otherwise).
+%
+% $Revision: 1.1.1.1 $ $Date: 2013/27/03 11:25:52 $
+% Original Author: Jorge Luis Bernal Rusiel
+% CVS Revision Info:
+% $Author: jbernal$
+% $Date: 2013/27/03 11:25:52 $
+% $Revision: 1.1 $
+% References: Kleinbaum, D.G., Klein, M., 2005. Survival analysis. A self-
+% learning approach, second edition. New York: Springer..
+%
+if nargin < 4
+ error('Too few inputs');
+elseif nargin < 5
+ e = 0.001;
+end;
+tic;
+[n,p] = size(X_ext);
+if (length(sID_ext)~=n) || (length(d_ext)~=n) || (length(t_ext)~=n)
+ error(['The design matrix X_ext, censorship status vector d_ext, time'...
+ ' vector t_ext and subject ID vector sID_ext must all have the same'...
+ ' number of rows.']);
+end;
+%indices of unique failure times in ft_ix (last index when ties happen)
+st_ix = find(d_ext==1);
+[~,ft_ix] = unique(t_ext(st_ix),'last');
+ft_ix = st_ix(ft_ix);
+%Starting values
+Bhat = zeros(p,1);
+
+%% Iterations
+nit = 50;
+gnorm = e+1;
+it = 1;
+display('Starting Newton-Raphson iterations');
+while (gnorm>e) && (it<=nit)
+ gr = SStat_Gradient(X_ext,t_ext,Bhat,ft_ix);
+ He = SStat_Hessian(X_ext,t_ext,Bhat,ft_ix);
+ if (cond(He) < 1e+10)
+ invHe = He\eye(p);
+ else
+ [Vtemp,Dtemp] = eig(He);
+ invHe = Vtemp*diag(1./max(diag(Dtemp),1e-5))*Vtemp';
+ end
+ Bhat = Bhat - invHe*gr;
+ %log-likelihood
+ llh = SStat_Likelihood(X_ext,t_ext,Bhat,ft_ix);
+ display(['Likelihood at iteration ' num2str(it) ' : ' num2str(llh)]);
+ gnorm = norm(gr);
+ display(['Gradient norm: ' num2str(gnorm)]);
+ it = it+1;
+end;
+%% Termination
+z_sc = Bhat./sqrt(diag(-invHe));
+pval = 2*(1-normcdf(abs(z_sc),0,1));
+stats = struct('Bhat',Bhat,'zscore',z_sc,'pval',pval,'CovBhat',-invHe,'llh',llh);
+if (gnorm<=e)
+ st = 1;
+else
+ st = 0;
+ display(['Algorithm does not converge after ' num2str(nit)...
+ ' iterations!!!']);
+end;
+et = toc;
+display(['Total elapsed time is ' num2str(et) ' seconds']);
+end
+
+
+
+
+
+
+%% Likelihood, Gradient and Hessian
+
+function llk = SStat_Likelihood(X_ext,t_ext,Bhat,ft_ix)
+%
+% Log-likelihood value.
+%
+% Input
+% X_ext: Extended design matrix.
+% t_ext: Extended survival time vector.
+% Bhat: Estimated vector of the population regression parameters.
+% ft_ix: Failure time indices in X_ext (last index if any tie).
+%
+% Output
+% llk: Log-likelihood value.
+%
+llk = 0;
+nft = length(ft_ix);
+for j=1:nft
+ term = 0;
+ ties = ft_ix(t_ext(ft_ix(j))==t_ext(ft_ix));
+ nties = length(ties)-1;
+ lpr = X_ext(ties,:)*Bhat;
+ aux1 = sum(exp(lpr));
+ aux2 = sum(exp(X_ext(t_ext(ft_ix(j))==t_ext,:)*Bhat));
+ for l=0:nties
+ term = term + log(aux2-l*aux1/(nties+1));
+ end;
+ llk = llk + sum(lpr) - term;
+end;
+end
+
+
+function gr = SStat_Gradient(X_ext,t_ext,Bhat,ft_ix)
+%
+% Gradient vector for the log-likelihood.
+%
+% Input
+% X_ext: Extended design matrix.
+% t_ext: Extended survival time vector.
+% Bhat: Estimated vector of the population regression parameters.
+% ft_ix: Failure time indices in X_ext (last index if any tie).
+%
+% Output
+% gr: Gradient vector.
+%
+p = size(X_ext,2);
+gr = zeros(p,1);
+nft = length(ft_ix);
+for j=1:nft
+ ties = ft_ix(t_ext(ft_ix(j))==t_ext(ft_ix));
+ nties = length(ties)-1;
+ riskset = t_ext(ft_ix(j))==t_ext;
+ term1 = sum(X_ext(ties,:),1);
+ term2 = exp(X_ext(riskset,:)*Bhat);
+ term3 = exp(X_ext(ties,:)*Bhat);
+ term4 = term2'*X_ext(riskset,:);
+ term5 = term3'*X_ext(ties,:);
+ term = 0;
+ for l=0:nties
+ term = term + (term4-l*term5/(nties+1))/(sum(term2)-l*sum(term3)/(nties+1));
+ end;
+ gr = gr + (term1 - term)';
+end;
+end
+
+
+function He = SStat_Hessian(X_ext,t_ext,Bhat,ft_ix)
+%
+% Hessian matrix for the log-likelihood.
+%
+% Input
+% X_ext: Extended design matrix.
+% t_ext: Extended survival time vector.
+% Bhat: Estimated vector of the population regression parameters.
+% ft_ix: Failure time indices in X_ext (last index if any tie).
+%
+% Output
+% He: Hessian matrix.
+%
+p = size(X_ext,2);
+He = zeros(p,p);
+nft = length(ft_ix);
+for j=1:nft
+ ties = ft_ix(t_ext(ft_ix(j))==t_ext(ft_ix));
+ nties = length(ties)-1;
+ rsk_ix = find(t_ext(ft_ix(j))==t_ext);
+ m = length(rsk_ix);
+ term1 = 0;
+ for i=1:m
+ term1 = term1 + exp(X_ext(rsk_ix(i),:)*Bhat)*X_ext(rsk_ix(i),:)'*X_ext(rsk_ix(i),:);
+ end;
+ term2 = 0;
+ for i=1:nties
+ term2 = term2 + exp(X_ext(ties(i),:)*Bhat)*X_ext(ties(i),:)'*X_ext(ties(i),:);
+ end;
+ term3 = sum(exp(X_ext(rsk_ix,:)*Bhat));
+ term4 = sum(exp(X_ext(ties,:)*Bhat));
+ term5 = exp(X_ext(rsk_ix,:)*Bhat)'*X_ext(rsk_ix,:);
+ term6 = exp(X_ext(ties,:)*Bhat)'*X_ext(ties,:);
+ term = 0;
+ for l=0:nties
+ Z = term5 - l*term6/(nties+1);
+ phi = term3 - l*term4/(nties+1);
+ term = term + (term1-l*term2/(nties+1))/phi - Z'*Z/(phi*phi);
+ end;
+ He = He - term;
+end;
+end
+
diff --git a/IQM Tools/IQMlite/tools/survival/IQMcoxPH.m b/IQM Tools/IQMlite/tools/survival/IQMcoxPH.m
new file mode 100644
index 0000000000000000000000000000000000000000..84cd86923744232a0917ceedd5b43cc1bc3f5913
--- /dev/null
+++ b/IQM Tools/IQMlite/tools/survival/IQMcoxPH.m
@@ -0,0 +1,241 @@
+function [stats,st] = IQMcoxPH(X,d,t,e)
+% [stats,st] = IQMcoxPH(X,d,t,e)
+%
+% Parameter estimation for the Cox proportional hazards model. Survival
+% time ties are handled using Efron's method.
+%
+% Input
+% X: Design Matrix with the time-independent covariates. (mxp, m # of
+% subjects, p # of covariates).
+% d: Logical vector (mx1) indicating censorship status (1 if the subject got
+% the failure event or 0 otherwise).
+% t: Vector (mx1) whose entries are the survival and censored times (ordered
+% according to X).
+% e: Convergence epsilon (gradient's norm). Default 10^-3;
+%
+% Output
+% stats.Bhat: Estimated vector of the population regression parameters.
+% stats.CovBhat: Estimated covariance matrix of the population regression
+% parameters.
+% stats.llh: Values of the maximum log-likelihood across the optimization
+% process.
+% st: Termination state (1 for convergence and 0 otherwise).
+%
+% $Revision: 1.1.1.1 $ $Date: 2013/27/03 11:25:52 $
+% Original Author: Jorge Luis Bernal Rusiel
+% CVS Revision Info:
+% $Author: jbernal$
+% $Date: 2013/27/03 11:25:52 $
+% $Revision: 1.1 $
+% References: Kleinbaum, D.G., Klein, M., 2005. Survival analysis. A self-
+% learning approach, second edition. New York: Springer..
+%
+if nargin < 3
+ error('Too few inputs');
+elseif nargin < 4
+ e = 0.001;
+end;
+tic;
+[m,p] = size(X);
+if (length(d)~=m) || (length(t)~=m)
+ error(['All, the design matrix X, censorship status vector d and'...
+ ' time vector t must have the same number of rows.']);
+end;
+%Sort the data by time. If there is a tie between a failure time and a
+%censored time then the failure time goes first.
+st_ix = find(d==1);
+t1 = t(st_ix);
+[t1,t1_ix] = sort(t1);
+X1 = X(st_ix(t1_ix),:);
+cs_ix = find(d==0);
+if ~isempty(cs_ix)
+ t2 = t(cs_ix);
+ [t2,t2_ix] = sort(t2);
+ X2 = X(cs_ix(t2_ix),:);
+ count1 = 1; count2 = 1; i = 0;
+ while (count1 <= length(t1)) && (count2 <= length(t2))
+ i = i + 1;
+ if t1(count1) <= t2(count2)
+ X(i,:) = X1(count1,:);
+ d(i) = 1;
+ t(i) = t1(count1);
+ count1 = count1 + 1;
+ else
+ X(i,:) = X2(count2,:);
+ d(i) = 0;
+ t(i) = t2(count2);
+ count2 = count2 + 1;
+ end;
+ end;
+ if (count1 > length(t1))
+ X(i+1:end,:) = X2(count2:end,:);
+ d(i+1:end) = 0;
+ t(i+1:end) = t2(count2:end);
+ else
+ X(i+1:end,:) = X1(count1:end,:);
+ d(i+1:end) = 1;
+ t(i+1:end) = t1(count1:end);
+ end;
+else
+ X = X1;
+ t = t1;
+end;
+%indices of unique failure times in ft_ix (last index when ties happen)
+st_ix = find(d==1);
+[ft,ft_ix] = unique(t(st_ix),'last');
+ft_ix = st_ix(ft_ix);
+nft = length(ft);
+%handling ties in failure times
+nties = zeros(1,nft);
+for j=1:nft
+ i = 1;
+ while (ft_ix(j)-i>0) && (ft(j)==t(ft_ix(j)-i))
+ i = i + 1;
+ end;
+ nties(j) = i-1;
+end;
+%Starting values
+Bhat = zeros(p,1);
+
+%% Iterations
+nit = 50;
+gnorm = e+1;
+it = 1;
+display('Starting Newton-Raphson iterations');
+while (gnorm>e) && (it<=nit)
+ gr = SStat_Gradient(X,Bhat,ft_ix,nties);
+ He = SStat_Hessian(X,Bhat,ft_ix,nties);
+ if (cond(He) < 1e+10)
+ invHe = He\eye(p);
+ else
+ [Vtemp,Dtemp] = eig(He);
+ invHe = Vtemp*diag(1./max(diag(Dtemp),1e-5))*Vtemp';
+ end
+ Bhat = Bhat - invHe*gr;
+ %log-likelihood
+ llh = SStat_Likelihood(X,Bhat,ft_ix,nties);
+ display(['Likelihood at iteration ' num2str(it) ' : ' num2str(llh)]);
+ gnorm = norm(gr);
+ display(['Gradient norm: ' num2str(gnorm)]);
+ it = it+1;
+end;
+%% Termination
+z_sc = Bhat./sqrt(diag(-invHe));
+pval = 2*(1-normcdf(abs(z_sc),0,1));
+stats = struct('Bhat',Bhat,'zscore',z_sc,'pval',pval,'CovBhat',-invHe,'llh',llh);
+if (gnorm<=e)
+ st = 1;
+else
+ st = 0;
+ display(['Algorithm does not converge after ' num2str(nit)...
+ ' iterations!!!']);
+end;
+et = toc;
+display(['Total elapsed time is ' num2str(et) ' seconds']);
+end
+
+
+
+
+
+
+%% Likelihood, Gradient and Hessian
+
+function llk = SStat_Likelihood(X,Bhat,ft_ix,nties)
+%
+% Log-likelihood value.
+%
+% Input
+% X: Ordered design Matrix (according to survival time).
+% Bhat: Estimated vector of the population regression parameters.
+% ft_ix: Failure time indices in X (last index if any tie).
+% nties: Number of ties for each survival time.
+%
+% Output
+% llk: Log-likelihood value.
+%
+llk = 0;
+nft = length(ft_ix);
+for j=1:nft
+ term = 0;
+ lpr = X(ft_ix(j)-nties(j):ft_ix(j),:)*Bhat;
+ aux1 = sum(exp(lpr));
+ aux2 = sum(exp(X(ft_ix(j)-nties(j):end,:)*Bhat));
+ for l=0:nties(j)
+ term = term + log(aux2-l*aux1/(nties(j)+1));
+ end;
+ llk = llk + sum(lpr) - term;
+end;
+end
+
+
+function gr = SStat_Gradient(X,Bhat,ft_ix,nties)
+%
+% Gradient vector for the log-likelihood.
+%
+% Input
+% X: Ordered design Matrix (according to survival time).
+% Bhat: Estimated vector of the population regression parameters.
+% ft_ix: Failure time indices in X (last index if any tie).
+% nties: Number of ties for each survival time.
+%
+% Output
+% gr: Gradient vector.
+%
+p = size(X,2);
+gr = zeros(p,1);
+nft = length(ft_ix);
+for j=1:nft
+ term1 = sum(X(ft_ix(j)-nties(j):ft_ix(j),:),1);
+ term2 = exp(X(ft_ix(j)-nties(j):end,:)*Bhat);
+ term3 = exp(X(ft_ix(j)-nties(j):ft_ix(j),:)*Bhat);
+ term4 = term2'*X(ft_ix(j)-nties(j):end,:);
+ term5 = term3'*X(ft_ix(j)-nties(j):ft_ix(j),:);
+ term = 0;
+ for l=0:nties(j)
+ term = term + (term4-l*term5/(nties(j)+1))/(sum(term2)-l*sum(term3)/(nties(j)+1));
+ end;
+ gr = gr + (term1 - term)';
+end;
+end
+
+
+function He = SStat_Hessian(X,Bhat,ft_ix,nties)
+%
+% Hessian matrix for the log-likelihood.
+%
+% Input
+% X: Ordered design Matrix (according to survival time).
+% Bhat: Estimated vector of the population regression parameters.
+% ft_ix: Failure time indices in X (last index if any tie).
+% nties: Number of ties for each survival time.
+%
+% Output
+% He: Hessian matrix.
+%
+[m,p] = size(X);
+He = zeros(p,p);
+nft = length(ft_ix);
+for j=1:nft
+ term1 = 0;
+ for i=ft_ix(j)-nties(j):m
+ term1 = term1 + exp(X(i,:)*Bhat)*X(i,:)'*X(i,:);
+ end;
+ term2 = 0;
+ for i=ft_ix(j)-nties(j):ft_ix(j)
+ term2 = term2 + exp(X(i,:)*Bhat)*X(i,:)'*X(i,:);
+ end;
+ term3 = sum(exp(X(ft_ix(j)-nties(j):end,:)*Bhat));
+ term4 = sum(exp(X(ft_ix(j)-nties(j):ft_ix(j),:)*Bhat));
+ term5 = exp(X(ft_ix(j)-nties(j):m,:)*Bhat)'*X(ft_ix(j)-nties(j):m,:);
+ term6 = exp(X(ft_ix(j)-nties(j):ft_ix(j),:)*Bhat)'*X(ft_ix(j)-nties(j):ft_ix(j),:);
+ term = 0;
+ for l=0:nties(j)
+ Z = term5 - l*term6/(nties(j)+1);
+ phi = term3 - l*term4/(nties(j)+1);
+ term = term + (term1-l*term2/(nties(j)+1))/phi - Z'*Z/(phi*phi);
+ end;
+ He = He - term;
+end;
+end
+
diff --git a/IQM Tools/IQMlite/tools/survival/IQMcoxStratPH.m b/IQM Tools/IQMlite/tools/survival/IQMcoxStratPH.m
new file mode 100644
index 0000000000000000000000000000000000000000..c63902d42e0f9ff3d612441e65ff075e8bbe8dd3
--- /dev/null
+++ b/IQM Tools/IQMlite/tools/survival/IQMcoxStratPH.m
@@ -0,0 +1,265 @@
+function [stats,st] = IQMcoxStratPH(X,z,d,t,e)
+% [stats,st] = IQMcoxPH(X,z,d,t,e)
+%
+% Parameter estimation for the Stratified Cox proportional hazards model.
+% Survival time ties are handled using Efron's method.
+%
+% Input
+% X: Design Matrix with the time-independent covariates. (mxp, m # of
+% subjects, p # of covariates).
+% z: Vector (mx1) with the strata for the rows of X;
+% d: Logical vector (mx1) indicating censorship status (1 if the subject got
+% the failure event or 0 otherwise).
+% t: Vector (mx1) whose entries are the survival and censored times (ordered
+% according to X).
+% e: Convergence epsilon (gradient's norm). Default 10^-3;
+%
+% Output
+% stats.Bhat: Estimated vector of the population regression parameters.
+% stats.CovBhat: Estimated covariance matrix of the population regression
+% parameters.
+% stats.llh: Values of the maximum log-likelihood across the optimization
+% process.
+% st: Termination state (1 for convergence and 0 otherwise).
+%
+% $Revision: 1.1.1.1 $ $Date: 2013/27/03 11:25:52 $
+% Original Author: Jorge Luis Bernal Rusiel
+% CVS Revision Info:
+% $Author: jbernal$
+% $Date: 2013/27/03 11:25:52 $
+% $Revision: 1.1 $
+% References: Kleinbaum, D.G., Klein, M., 2005. Survival analysis. A self-
+% learning approach, second edition. New York: Springer..
+%
+if nargin < 4
+ error('Too few inputs');
+elseif nargin < 5
+ e = 0.001;
+end;
+tic;
+[m,p] = size(X);
+if (length(z)~=m) || (length(d)~=m) || (length(t)~=m)
+ error(['The design matrix X, strata vector z, censorship status vector d'...
+ ' and time vector t must all have the same number of rows.']);
+end;
+%Sort the data by time. If there is a tie between a failure time and a
+%censored time then the failure time goes first.
+st_ix = find(d==1);
+t1 = t(st_ix);
+[t1,t1_ix] = sort(t1);
+X1 = X(st_ix(t1_ix),:);
+z1 = z(st_ix(t1_ix));
+cs_ix = find(d==0);
+if ~isempty(cs_ix)
+ t2 = t(cs_ix);
+ [t2,t2_ix] = sort(t2);
+ X2 = X(cs_ix(t2_ix),:);
+ z2 = z(cs_ix(t2_ix));
+ count1 = 1; count2 = 1; i = 0;
+ while (count1 <= length(t1)) && (count2 <= length(t2))
+ i = i + 1;
+ if t1(count1) <= t2(count2)
+ X(i,:) = X1(count1,:);
+ z(i) = z1(count1);
+ d(i) = 1;
+ t(i) = t1(count1);
+ count1 = count1 + 1;
+ else
+ X(i,:) = X2(count2,:);
+ z(i) = z2(count2);
+ d(i) = 0;
+ t(i) = t2(count2);
+ count2 = count2 + 1;
+ end;
+ end;
+ if (count1 > length(t1))
+ X(i+1:end,:) = X2(count2:end,:);
+ z(i+1:end) = z2(count2:end);
+ d(i+1:end) = 0;
+ t(i+1:end) = t2(count2:end);
+ else
+ X(i+1:end,:) = X1(count1:end,:);
+ z(i+1:end) = z1(count1:end);
+ d(i+1:end) = 1;
+ t(i+1:end) = t1(count1:end);
+ end;
+else
+ X = X1;
+ t = t1;
+ z = z1;
+end;
+uz = unique(z);
+nst = length(uz);
+%indices of unique failure times for each stratum in ft_ix (last index when ties happen)
+ft_ix = zeros(nst,sum(d));
+ft = ft_ix;
+nft = zeros(nst,1);
+nties = ft;
+for k=1:nst
+ t_strat = t(z==uz(k));
+ st_ix = find(d(z==uz(k))==1);
+ [ft_aux,ft_ix_aux] = unique(t_strat(st_ix),'last');
+ nft(k) = length(ft_aux);
+ ft_ix(k,1:nft(k)) = st_ix(ft_ix_aux);
+ ft(k,1:nft(k)) = ft_aux;
+ %handling ties in failure times
+ for j=1:nft(k)
+ i = 1;
+ while (ft_ix(k,j)-i>0) && (ft(k,j)==t(ft_ix(k,j)-i))
+ i = i + 1;
+ end;
+ nties(k,j) = i-1;
+ end;
+end;
+%Starting values
+Bhat = zeros(p,1);
+
+%% Iterations
+nit = 50;
+gnorm = e+1;
+it = 1;
+display('Starting Newton-Raphson iterations');
+while (gnorm>e) && (it<=nit)
+ gr = zeros(p,1);
+ He = zeros(p,p);
+ for k=1:length(uz)
+ gr = gr + SStat_Gradient(X(z==uz(k),:),Bhat,ft_ix(k,1:nft(k)),nties(k,1:nft(k)));
+ He = He + SStat_Hessian(X(z==uz(k),:),Bhat,ft_ix(k,1:nft(k)),nties(k,1:nft(k)));
+ end;
+ if (cond(He) < 1e+10)
+ invHe = He\eye(p);
+ else
+ [Vtemp,Dtemp] = eig(He);
+ invHe = Vtemp*diag(1./max(diag(Dtemp),1e-5))*Vtemp';
+ end
+ Bhat = Bhat - invHe*gr;
+ %log-likelihood
+ llh = 0;
+ for k=1:length(uz)
+ llh = llh + SStat_Likelihood(X(z==uz(k),:),Bhat,ft_ix(k,1:nft(k)),nties(k,1:nft(k)));
+ end;
+ display(['Likelihood at iteration ' num2str(it) ' : ' num2str(llh)]);
+ gnorm = norm(gr);
+ display(['Gradient norm: ' num2str(gnorm)]);
+ it = it+1;
+end;
+%% Termination
+z_sc = Bhat./sqrt(diag(-invHe));
+pval = 2*(1-normcdf(abs(z_sc),0,1));
+stats = struct('Bhat',Bhat,'zscore',z_sc,'pval',pval,'CovBhat',-invHe,'llh',llh);
+if (gnorm<=e)
+ st = 1;
+else
+ st = 0;
+ display(['Algorithm does not converge after ' num2str(nit)...
+ ' iterations!!!']);
+end;
+et = toc;
+display(['Total elapsed time is ' num2str(et) ' seconds']);
+end
+
+
+
+
+
+
+%% Likelihood, Gradient and Hessian
+
+function llk = SStat_Likelihood(X,Bhat,ft_ix,nties)
+%
+% Log-likelihood value.
+%
+% Input
+% X: Ordered design Matrix (according to survival time).
+% Bhat: Estimated vector of the population regression parameters.
+% ft_ix: Failure time indices in X (last index if any tie).
+% nties: Number of ties for each survival time.
+%
+% Output
+% llk: Log-likelihood value.
+%
+llk = 0;
+nft = length(ft_ix);
+for j=1:nft
+ term = 0;
+ lpr = X(ft_ix(j)-nties(j):ft_ix(j),:)*Bhat;
+ aux1 = sum(exp(lpr));
+ aux2 = sum(exp(X(ft_ix(j)-nties(j):end,:)*Bhat));
+ for l=0:nties(j)
+ term = term + log(aux2-l*aux1/(nties(j)+1));
+ end;
+ llk = llk + sum(lpr) - term;
+end;
+end
+
+
+function gr = SStat_Gradient(X,Bhat,ft_ix,nties)
+%
+% Gradient vector for the log-likelihood.
+%
+% Input
+% X: Ordered design Matrix (according to survival time).
+% Bhat: Estimated vector of the population regression parameters.
+% ft_ix: Failure time indices in X (last index if any tie).
+% nties: Number of ties for each survival time.
+%
+% Output
+% gr: Gradient vector.
+%
+p = size(X,2);
+gr = zeros(p,1);
+nft = length(ft_ix);
+for j=1:nft
+ term1 = sum(X(ft_ix(j)-nties(j):ft_ix(j),:),1);
+ term2 = exp(X(ft_ix(j)-nties(j):end,:)*Bhat);
+ term3 = exp(X(ft_ix(j)-nties(j):ft_ix(j),:)*Bhat);
+ term4 = term2'*X(ft_ix(j)-nties(j):end,:);
+ term5 = term3'*X(ft_ix(j)-nties(j):ft_ix(j),:);
+ term = 0;
+ for l=0:nties(j)
+ term = term + (term4-l*term5/(nties(j)+1))/(sum(term2)-l*sum(term3)/(nties(j)+1));
+ end;
+ gr = gr + (term1 - term)';
+end;
+end
+
+
+function He = SStat_Hessian(X,Bhat,ft_ix,nties)
+%
+% Hessian matrix for the log-likelihood.
+%
+% Input
+% X: Ordered design Matrix (according to survival time).
+% Bhat: Estimated vector of the population regression parameters.
+% ft_ix: Failure time indices in X (last index if any tie).
+% nties: Number of ties for each survival time.
+%
+% Output
+% He: Hessian matrix.
+%
+[m,p] = size(X);
+He = zeros(p,p);
+nft = length(ft_ix);
+for j=1:nft
+ term1 = 0;
+ for i=ft_ix(j)-nties(j):m
+ term1 = term1 + exp(X(i,:)*Bhat)*X(i,:)'*X(i,:);
+ end;
+ term2 = 0;
+ for i=ft_ix(j)-nties(j):ft_ix(j)
+ term2 = term2 + exp(X(i,:)*Bhat)*X(i,:)'*X(i,:);
+ end;
+ term3 = sum(exp(X(ft_ix(j)-nties(j):end,:)*Bhat));
+ term4 = sum(exp(X(ft_ix(j)-nties(j):ft_ix(j),:)*Bhat));
+ term5 = exp(X(ft_ix(j)-nties(j):m,:)*Bhat)'*X(ft_ix(j)-nties(j):m,:);
+ term6 = exp(X(ft_ix(j)-nties(j):ft_ix(j),:)*Bhat)'*X(ft_ix(j)-nties(j):ft_ix(j),:);
+ term = 0;
+ for l=0:nties(j)
+ Z = term5 - l*term6/(nties(j)+1);
+ phi = term3 - l*term4/(nties(j)+1);
+ term = term + (term1-l*term2/(nties(j)+1))/phi - Z'*Z/(phi*phi);
+ end;
+ He = He - term;
+end;
+end
+
diff --git a/IQM Tools/IQMlite/tools/survival/IQMhr.m b/IQM Tools/IQMlite/tools/survival/IQMhr.m
new file mode 100644
index 0000000000000000000000000000000000000000..018d5d30adf7fd2496fba198dceb471566e8dadc
--- /dev/null
+++ b/IQM Tools/IQMlite/tools/survival/IQMhr.m
@@ -0,0 +1,41 @@
+function [hr,pval,CI] = IQMhr(x1,x2,stats)
+% [hr,pval,CI] = IQMhr(x1,x2,stats)
+%
+% Hazard ratio estimate for Cox models.
+%
+% Input
+%
+% x1: Row vector with the covariate values for the first group.
+% x2: Row vector with the covariate values for the second group.
+% stats: Structure containing statistiscs obtained with any of IQMcoxPH,
+% IQMcoxStratPH and IQMcoxExt.
+%
+% Output
+% hr: Hazard ratio value.
+% pval: P-value for the hr value.
+% CI: 95% confidence interval for the hazard ratio estimate.
+%
+% $Revision: 1.1.1.1 $ $Date: 2013/27/03 11:25:52 $
+% Original Author: Jorge Luis Bernal Rusiel
+% CVS Revision Info:
+% $Author: jbernal$
+% $Date: 2013/27/03 11:25:52 $
+% $Revision: 1.1 $
+% References: Kleinbaum, D.G., Klein, M., 2005. Survival analysis. A self-
+% learning approach, second edition. New York: Springer..
+%
+if nargin < 3
+ error('Too few inputs');
+end;
+p = length(stats.Bhat);
+if (length(x1)~=p) || (length(x2)~=p)
+ error(['Vectors x1, x2 and stats.Bhat must have the same length.']);
+end;
+d = x1-x2;
+lrp = d*stats.Bhat;
+hr = exp(lrp);
+lrp_std = (d'*d).*stats.CovBhat;
+lrp_std = sqrt(sum(lrp_std(:)));
+CI.low = exp(lrp-1.96*lrp_std);
+CI.high = exp(lrp+1.96*lrp_std);
+pval = 2*(1-normcdf(abs(lrp/lrp_std),0,1));
\ No newline at end of file
diff --git a/IQM Tools/IQMlite/tools/survival/IQMkm.m b/IQM Tools/IQMlite/tools/survival/IQMkm.m
new file mode 100644
index 0000000000000000000000000000000000000000..9a0e5b80121ab56a02e316a6b6adc4144f103c9c
--- /dev/null
+++ b/IQM Tools/IQMlite/tools/survival/IQMkm.m
@@ -0,0 +1,59 @@
+function [km] = IQMkm(z,d,t,ptype)
+% [km] = IQMkm(z,d,t,ptype)
+%
+% Kaplan-Meier curves estimation and plot.
+%
+% Input
+% z: Discrete covariate (mx1, m # of subjects) with the strata categories.
+% d: Logical vector (mx1) indicating censorship status (1 if the subject got
+% the failure event or 0 otherwise).
+% t: Vector (mx1) whose entries are the survival and censored times (ordered
+% according to X).
+% ptype: Plot type string. It can only be either 'km' to plot the KM curves
+% or 'll' to plot the -log(-log(KM)) curves.
+%
+% Output
+% km: Structure array with the Kaplan-Meier estimates for each category of
+% z.
+%
+% $Revision: 1.1.1.1 $ $Date: 2013/27/03 11:25:52 $
+% Original Author: Jorge Luis Bernal Rusiel
+% CVS Revision Info:
+% $Author: jbernal$
+% $Date: 2013/27/03 11:25:52 $
+% $Revision: 1.1 $
+% References: Kleinbaum, D.G., Klein, M., 2005. Survival analysis. A self-
+% learning approach, second edition. New York: Springer.
+%
+if nargin < 3
+ error('Too few inputs');
+elseif nargin < 4
+ ptype = 'km';
+end;
+uz = unique(z);
+nuz = length(uz);
+km(1:nuz) = struct('f',[],'t',[]);
+if strcmpi(ptype,'km')
+ figure('Name','KM plot');
+elseif strcmpi(ptype,'ll')
+ figure('Name','-log(-log(KM)) plot');
+else
+ error('The only posible values for ptype are ''km'' or ''ll''.');
+end;
+hold on
+colors = {'b','k','g','c','m','r','y'};
+for i=1:length(uz)
+ t_strat = t(z==uz(i));
+ d_strat = d(z==uz(i));
+ [km(i).f,km(i).t] = ecdf(t_strat,'censoring',~d_strat);
+ km(i).f = 1 - km(i).f;
+ if strcmpi(ptype,'km')
+ stairs(km(i).t,km(i).f,colors{mod(i,7)+1},'LineWidth',2);
+ else
+ stairs(km(i).t,-log(-log(km(i).f)),colors{mod(i,7)+1},'LineWidth',2);
+ end;
+end;
+legend(num2str(uz));
+hold off
+
+
diff --git a/IQM Tools/IQMlite/tools/survival/IQMphTest.m b/IQM Tools/IQMlite/tools/survival/IQMphTest.m
new file mode 100644
index 0000000000000000000000000000000000000000..2c427c49f975347b3368e4266274a5e3df2fb3e3
--- /dev/null
+++ b/IQM Tools/IQMlite/tools/survival/IQMphTest.m
@@ -0,0 +1,131 @@
+function [pval,rho] = IQMphTest(X,d,t,stats_PH)
+% [pval,rho] = IQMphTest(X,d,t,stats_PH)
+%
+% Schoenfeld residuals test for the proportional hazards assumption.
+%
+% Input
+% X: Design Matrix with the time-independent covariates. (mxp, m # of
+% subjects, p # of covariates).
+% d: Logical vector (mx1) indicating censorship status (1 if the subject got
+% the failure event or 0 otherwise).
+% t: Vector (mx1) whose entries are the survival and censored times (ordered
+% according to X).
+% stats_PH: Structure obtained from IQMcoxPH
+%
+% Output
+% pval: P-values indicating the probability of the PH assumption for each
+% covariate.
+% rho: Correlations between Schoenfeld residuals and ranked failure times.
+%
+% $Revision: 1.1.1.1 $ $Date: 2013/27/03 11:25:52 $
+% Original Author: Jorge Luis Bernal Rusiel
+% CVS Revision Info:
+% $Author: jbernal$
+% $Date: 2013/27/03 11:25:52 $
+% $Revision: 1.1 $
+% References: Kleinbaum, D.G., Klein, M., 2005. Survival analysis. A self-
+% learning approach, second edition. New York: Springer..
+%
+if nargin < 4
+ error('Too few inputs');
+end;
+[m,p] = size(X);
+if (length(d)~=m) || (length(t)~=m)
+ error(['All, the design matrix X, censorship status vector d and'...
+ ' time vector t must have the same number of rows.']);
+end;
+if length(stats_PH.Bhat)~=p
+ error(['The vector stats.Bhat must have the same length as the number'...
+ ' of colums of X.']);
+end;
+%Sort the data by time. If there is a tie between a failure time and a
+%censored time then the failure time goes first.
+st_ix = find(d==1);
+t1 = t(st_ix);
+[t1,t1_ix] = sort(t1);
+X1 = X(st_ix(t1_ix),:);
+cs_ix = find(d==0);
+if ~isempty(cs_ix)
+ t2 = t(cs_ix);
+ [t2,t2_ix] = sort(t2);
+ X2 = X(cs_ix(t2_ix),:);
+ count1 = 1; count2 = 1; i = 0;
+ while (count1 <= length(t1)) && (count2 <= length(t2))
+ i = i + 1;
+ if t1(count1) <= t2(count2)
+ X(i,:) = X1(count1,:);
+ d(i) = 1;
+ t(i) = t1(count1);
+ count1 = count1 + 1;
+ else
+ X(i,:) = X2(count2,:);
+ d(i) = 0;
+ t(i) = t2(count2);
+ count2 = count2 + 1;
+ end;
+ end;
+ if (count1 > length(t1))
+ X(i+1:end,:) = X2(count2:end,:);
+ d(i+1:end) = 0;
+ t(i+1:end) = t2(count2:end);
+ else
+ X(i+1:end,:) = X1(count1:end,:);
+ d(i+1:end) = 1;
+ t(i+1:end) = t1(count1:end);
+ end;
+else
+ X = X1;
+ t = t1;
+end;
+%indices of unique failure times in ft_ix (last index when ties happen)
+st_ix = find(d==1);
+[ft,ft_ix] = unique(t(st_ix),'last');
+ft_ix = st_ix(ft_ix);
+nft = length(ft);
+%handling ties in failure times by substracting a very small random number
+rand('state',sum(100*clock));
+for j=1:nft
+ i = 1;
+ while (ft_ix(j)-i>0) && (ft(j)==t(ft_ix(j)-i))
+ i = i + 1;
+ end;
+ nties = i-1;
+ tt = t(ft_ix(j)-nties:ft_ix(j)) - 1e-5*(1 + rand(nties+1,1));
+ tX = X(ft_ix(j)-nties:ft_ix(j),:);
+ [stt,stt_ix] = sort(tt);
+ t(ft_ix(j)-nties:ft_ix(j)) = stt;
+ X(ft_ix(j)-nties:ft_ix(j),:) = tX(stt_ix,:);
+end;
+nft = length(st_ix);
+r = zeros(nft,p);
+for i=1:nft
+ lprv = exp(X(st_ix(i):end,:)*stats_PH.Bhat);
+ r(i,:) = X(st_ix(i),:) - (lprv'*X(st_ix(i):end,:))./sum(lprv);
+end;
+[rho,pval] = corr(r,[1:nft]');
+
+% Plot Schoenfeld residuals over Ranked failure times for eye inspection.
+figure;
+nall = size(r,2);
+ncol = ceil(sqrt(nall));
+nrow = ceil(nall/ncol);
+for k=1:nall,
+ subplot(nrow,ncol,k);
+ plot([1:nft],r(:,k),'o');
+ ylabel('Schoenfeld residuals');
+ xlabel('Ranked failure times');
+ title(sprintf('beta_%d',k),'Interpreter','none');
+ grid on
+ set(gca,'XLim',[1 nft])
+ % Add loess line
+ hold on
+ rsmooth = smoothIQM([1:nft],r(:,k),0.5,'loess');
+ plot([1:nft],rsmooth,'k-','LineWidth',2);
+ % Plot zero line
+ plot(get(gca,'XLim'),[0 0],'k--');
+ % Add legend
+ h = legend(['Residuals over time ' sprintf('\n(corr=%1.2g, p=%1.2g)',rho(k),pval(k))],'Loess smoother','Zero line','Location','best');
+ set(h,'Interpreter','none')
+ axis square
+end
+
diff --git a/IQM Tools/IQMlite/tools/survival/IQMxext.m b/IQM Tools/IQMlite/tools/survival/IQMxext.m
new file mode 100644
index 0000000000000000000000000000000000000000..ae3fa25311d2160cd87b5b5a02d5ba26461e2e9a
--- /dev/null
+++ b/IQM Tools/IQMlite/tools/survival/IQMxext.m
@@ -0,0 +1,114 @@
+function [sID_ext,X_ext,d_ext,t_ext] = IQMxext(sID,X,d,t)
+% [sID_ext,X_ext,d_ext,t_ext] = IQMxext(sID,X,d,t)
+%
+% Extension (by adding new rows) of a design matrix comprising only time-
+% independent covariates to a design matrix for which time-dependent
+% covariates can be easily added.
+%
+% Input
+% sID: Subjects' IDs (for each row of X).
+% X: Design Matrix with the time-independent covariates. (mxp, m # of
+% subjects, p # of covariates).
+% d: Logical vector (mx1) indicating censorship status (1 if the subject got
+% the failure event or 0 otherwise).
+% t: Vector (mx1) whose entries are the survival and censored times (ordered
+% according to X).
+%
+% Output
+% sID_ext: Extended subjects' IDs.
+% X_ext: Extended design matrix.
+% d_ext: Extended censorship status vector.
+% t_ext: Extended survival time vector.
+%
+% $Revision: 1.1.1.1 $ $Date: 2013/27/03 11:25:52 $
+% Original Author: Jorge Luis Bernal Rusiel
+% CVS Revision Info:
+% $Author: jbernal$
+% $Date: 2013/27/03 11:25:52 $
+% $Revision: 1.1 $
+% References: Kleinbaum, D.G., Klein, M., 2005. Survival analysis. A self-
+% learning approach, second edition. New York: Springer..
+%
+if nargin < 4
+ error('Too few inputs');
+end;
+[m,p] = size(X);
+if (length(sID)~=m) || (length(d)~=m) || (length(t)~=m)
+ error(['All, the design matrix X, censorship status vector d, '...
+ 'time vector t and subject ID vector must have the same number of rows.']);
+end;
+%Sort the data by time. If there is a tie between a failure time and a
+%censored time then the failure time goes first.
+st_ix = find(d==1);
+t1 = t(st_ix);
+[t1,t1_ix] = sort(t1);
+X1 = X(st_ix(t1_ix),:);
+sID1 = sID(st_ix(t1_ix));
+cs_ix = find(d==0);
+if ~isempty(cs_ix)
+ t2 = t(cs_ix);
+ [t2,t2_ix] = sort(t2);
+ X2 = X(cs_ix(t2_ix),:);
+ sID2 = sID(cs_ix(t2_ix));
+ count1 = 1; count2 = 1; i = 0;
+ while (count1 <= length(t1)) && (count2 <= length(t2))
+ i = i + 1;
+ if t1(count1) <= t2(count2)
+ sID(i) = sID1(count1);
+ X(i,:) = X1(count1,:);
+ d(i) = 1;
+ t(i) = t1(count1);
+ count1 = count1 + 1;
+ else
+ sID(i) = sID2(count2);
+ X(i,:) = X2(count2,:);
+ d(i) = 0;
+ t(i) = t2(count2);
+ count2 = count2 + 1;
+ end;
+ end;
+ if (count1 > length(t1))
+ sID(i+1:end) = sID2(count2:end);
+ X(i+1:end,:) = X2(count2:end,:);
+ d(i+1:end) = 0;
+ t(i+1:end) = t2(count2:end);
+ else
+ sID(i+1:end) = sID1(count1:end);
+ X(i+1:end,:) = X1(count1:end,:);
+ d(i+1:end) = 1;
+ t(i+1:end) = t1(count1:end);
+ end;
+else
+ sID = sID1;
+ X = X1;
+ t = t1;
+end;
+%indices of unique failure times in ft_ix (last index when ties happen)
+st_ix = find(d==1);
+[ft,ft_ix] = unique(t(st_ix),'last');
+ft_ix = st_ix(ft_ix);
+nft = length(ft);
+n = 0;
+for j=1:nft
+ n = n + sum(t>t(ft_ix(j)));
+end;
+n = n + m;
+X_ext = zeros(n,p);
+d_ext = zeros(n,1);
+t_ext = zeros(n,1);
+count = 1;
+for i=1:m
+ ni = sum(t(ft_ix) < t(i));
+ X_ext(count:count+ni,:) = kron(ones(ni+1,1),X(i,:));
+ sID_ext(count:count+ni) = sID(i);
+ d_ext(count:count+ni-1) = 0;
+ d_ext(count+ni) = d(i);
+ t_ext(count:count+ni-1) = t(ft_ix(t(ft_ix)<t(i)));
+ t_ext(count+ni) = t(i);
+ count = count + ni + 1;
+end;
+
+
+
+
+
diff --git a/IQM Tools/IQMlite/tools/templates/scriptTemplateIQM.m b/IQM Tools/IQMlite/tools/templates/scriptTemplateIQM.m
index 4aec74ecd49ff198381bf1a05081d9669638186a..eedba36c5bce258c5779a6f805f7744f3935ec7d 100644
--- a/IQM Tools/IQMlite/tools/templates/scriptTemplateIQM.m
+++ b/IQM Tools/IQMlite/tools/templates/scriptTemplateIQM.m
@@ -1,43 +1,45 @@
-%% In this section the purpose of the script should be described.
-% Purpose of activity:
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% SCRIPT_NAME
%
-% Required previous scripts to be run:
+% [PURPOSE]
+%
+%
+% [AUTHOR]
%
-% Author:
%
-% Date:
+% [DATE]
+%
%
+% [DEPENDENCIES]
+%
+%
+% [INPUT]
+%
+%
+%
+% [OUTPUT]
+%
+% [OTHER]
+%
+%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-%% Installation of IQM Tools
-% This block should always be on the top of every analysis script with IQM
-% Tools. It might only be ommitted in the case that a central installation
-% of IQM Tools is present on the computer system and the IQM Tools are
-% automatically loaded during the startup of MATLAB.
-clc; % Clear command window
-clear all; % Clear workspace from all defined variables
-close all; % Close all figures
-fclose all; % Close all open files
-restoredefaultpath(); % Clear all user installed toolboxes
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% Clear all - install needed scripts
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+clc;
+clear all;
+close all;
+fclose all;
+restoredefaultpath;
% In the next line, please enter the path to your IQM Tools folder. It can
% be a relative or an absolute path.
PATH_IQM = 'D:\IQM Tools Suite';
-oldpath = pwd(); cd(PATH_IQM); installIQMtoolsInitial; cd(oldpath);
-
-%% "Initialize the compliance mode".
-% IQM Tools allows to automatically generate logfiles for all output that
-% is generated using the functions "IQMprintFigure" and
-% "IQMwriteText2File". These logfiles contain the username, the name of the
-% generated file (including the path), and the scripts and functions that
-% have been called to generate the output file. In order to ensure this is
-% working correctly, the only thing that needs to be done is to execute the
-% following command at the start of each analysis script.
-
-% The input argument to the "IQMinitializeCompliance" function needs to be
-% the name of the script file.
-IQMinitializeCompliance('TEMPLATE_SCRIPT_NAME');
-
-%% Here your code starts
-
+oldpath = pwd(); cd(PATH_IQM); installIQMtools; cd(oldpath);
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% Intialize compliance mode
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+IQMinitializeCompliance('SCRIPT_NAME.m');
diff --git a/IQM Tools/IQMpro/Contents.m b/IQM Tools/IQMpro/Contents.m
index 25a52d1094305764153b1bb5d2d3127a781fdd8c..c013ffb8abc60584abdb25a88e12e1e7e64699c4 100644
--- a/IQM Tools/IQMpro/Contents.m
+++ b/IQM Tools/IQMpro/Contents.m
@@ -1,5 +1,5 @@
% IQM Tools Pro
-% Version 1.2.1 (R2015B) 30.04.2016
+% Version 1.2.2 (R2015B) 02.01.2017
%
% IQMpro
% =======
diff --git a/IQM Tools/IQMpro/SETUP_PATHS_TOOLS_IQMPRO.m b/IQM Tools/IQMpro/SETUP_PATHS_TOOLS_IQMPRO.m
index 8cb89316d93e2e4a460b1ccd299a21015f101624..4fed82542cdade9bb73a9ac08fe1fbbb4c919d7b 100644
--- a/IQM Tools/IQMpro/SETUP_PATHS_TOOLS_IQMPRO.m
+++ b/IQM Tools/IQMpro/SETUP_PATHS_TOOLS_IQMPRO.m
@@ -13,16 +13,17 @@
% assumed that the command that runs a NONMEM or MONOLIX job via the queue
% returns the jobID as a number. It is also assumed that this number
% appears when the queue status command is called:
-PATH_SYSTEM_QUEUE_STATUS_UNIX = 'squeue';
+PATH_SYSTEM_QUEUE_STATUS_UNIX = '';
% If no queuing system used then keep this variable empty ('')
% NONMEM (currently tested version: 7.2 and 7.3)
PATH_SYSTEM_NONMEM_WINDOWS = 'nmfe73';
-PATH_SYSTEM_NONMEM_UNIX = 'qnmfe -n -t 1440 -x'; % Define UNIX path for Mac
+PATH_SYSTEM_NONMEM_UNIX = 'nmfe73';
% NONMEM PARALLEL
PATH_SYSTEM_NONMEM_PARALLEL_WINDOWS = 'nmfe73par';
-PATH_SYSTEM_NONMEM_PARALLEL_UNIX = 'qnmfe -n -t 1440 -x -c'; % Define UNIX path for Mac
+PATH_SYSTEM_NONMEM_PARALLEL_UNIX = 'nmfe73par';
+
% nmfe73par is assumed to be a shell script with the following calling
% syntax: "nmfe73par NRNODES controlfile outputfile". If you do not have
% one available, please ask your sysadmin to generate on for you
diff --git a/IQM Tools/IQMpro/auxiliary/data/checkDataColumnsIQM.m b/IQM Tools/IQMpro/auxiliary/data/checkDataColumnsIQM.m
index c973bef47126e1dc7833daee16a27f5185a53718..71b739ac696bb4c98a095b32749ba7f62981eef0 100644
--- a/IQM Tools/IQMpro/auxiliary/data/checkDataColumnsIQM.m
+++ b/IQM Tools/IQMpro/auxiliary/data/checkDataColumnsIQM.m
@@ -33,8 +33,10 @@ end
colnamestable = data.Properties.VariableNames;
errorText = '';
for k=1:length(COLNAMES),
- if isempty(strmatchIQM(COLNAMES{k},colnamestable,'exact')),
- errorText = sprintf('%sThe dataset does not contain the column "%s".\n',errorText,COLNAMES{k});
+ if ~isempty(COLNAMES{k}),
+ if isempty(strmatchIQM(COLNAMES{k},colnamestable,'exact')),
+ errorText = sprintf('%sThe dataset does not contain the column "%s".\n',errorText,COLNAMES{k});
+ end
end
end
diff --git a/IQM Tools/IQMpro/release.txt b/IQM Tools/IQMpro/release.txt
index 42c36a424641ae60bfd590acee81242f71a16029..7ab71b93e7796aa74aba9450f0a7ae6221b65a30 100644
--- a/IQM Tools/IQMpro/release.txt
+++ b/IQM Tools/IQMpro/release.txt
@@ -1,6 +1,11 @@
Release Notes IQM Tools Pro
===========================
+Version 1.2.2 (02.01.2017)
+--------------------------
+- Minor big fixes
+- Improved graphical output
+
Version 1.2.1 (30.04.2016)
--------------------------
- Minor bug fixes
diff --git a/IQM Tools/IQMpro/tools/01-MEXmodels/CVODEMEX/lib/CVODEmex25.a b/IQM Tools/IQMpro/tools/01-MEXmodels/CVODEMEX/lib/CVODEmex25.a
new file mode 100644
index 0000000000000000000000000000000000000000..fd4d2f3c6a2cac80df2cee1ca4437312bc150e7a
Binary files /dev/null and b/IQM Tools/IQMpro/tools/01-MEXmodels/CVODEMEX/lib/CVODEmex25.a differ
diff --git a/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMaddNLMEinfo2data.m b/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMaddNLMEinfo2data.m
index 5a9fc43af8d073ae801079f24066f4f3e1f0d70a..c82ed919110cc7578f30579a7f19e2abc1bff5a7 100644
--- a/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMaddNLMEinfo2data.m
+++ b/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMaddNLMEinfo2data.m
@@ -6,7 +6,10 @@ function [dataOut] = IQMaddNLMEinfo2data(data,DOSENAMES,OBSNAMES,FLAG_EXPAND_REP
%
% Repeated doses, defined by entries in columns INTERVAL and NRDOSES are
% NOT expanded by default. Expansion can be done by setting input argument
-% FLAG_EXPAND_REPEATED_DOSES = 1 (default: 0).
+% FLAG_EXPAND_REPEATED_DOSES = 1 (default: 0). Note that NRDOSES codes for
+% ADDITIONAL doses, in the same way (identical) as the ADDL column in
+% NONMEM and MONOLIX. This means that NRDOSES=1 codes for 2 doses, spaced
+% by INTERVAL.
%
% The following columns are added (if already present they are overwritten):
%
@@ -93,6 +96,9 @@ function [dataOut] = IQMaddNLMEinfo2data(data,DOSENAMES,OBSNAMES,FLAG_EXPAND_REP
% INTERVAL and NRDOSES are NOT expanded by default.
% Expansion can be done by setting input argument
% FLAG_EXPAND_REPEATED_DOSES = 1 (default: 0).
+% NOTE THAT NRDOSES has the same meaning as ADDL in
+% NONMEM and MONOLIX ... so it is the number of
+% ADDITIONAL doses.
%
% [OUTPUT]
% dataOut: Dataset with added NLME information.
diff --git a/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMcheckGeneralDataFormat.m b/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMcheckGeneralDataFormat.m
index 9ed29d754e6047ce3ac235ffaada74ff470aca0e..bf040aee8685286fa4caba4f8a2e4d49f37d55da 100644
--- a/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMcheckGeneralDataFormat.m
+++ b/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMcheckGeneralDataFormat.m
@@ -61,8 +61,11 @@ function [data] = IQMcheckGeneralDataFormat(data,silent)
% or the value of other readouts. The values need to be in the units, defined in the UNIT column.
% Specific cases:
% - For concomitant medications the dose will be given
-% - For BLOQ values, 0 will be used
% - Severity levels for adverse events
+% - For BLOQ records: any value can be entered that is lower than the actual LLOQ.
+% It is not acceptable to set this value to �NaN� or �NA� since then no discrimination
+% can be made between �missing� and �<LLOQ�. For PK records on untransformed data �0�
+% is suggested. On log transformed data 0 should not be used but log(LLOQ/2) would be acceptable.
% Should not be populated if VALUETXT is populated
% VALUETXT string Text version of value (if available and useful)
% Character value as given in the CRF.
@@ -72,7 +75,9 @@ function [data] = IQMcheckGeneralDataFormat(data,silent)
% LLOQ numeric Lower limit of quantification of event defined by NAME
% ROUTE string Route of administration (iv, subcut, intramuscular, intraarticular, oral, inhaled, topical, rectal)
% INTERVAL numeric Interval of dosing, if single row should define multiple dosings
-% NRDOSES numeric Number of doses given with the specified interval, if single row should define multiple dosings
+% NRDOSES numeric Number of ADDITIONAL doses given within the specified interval,
+% NRDOSES=N codes for a total of N+1 doses with an interval as defined in the
+% column "INTERVAL" starting at the time defined in the dose record.
% COMMENT string Additional information for the observation/event
% For example:
% - For PK: concatenation of DMPK flag to exclude or not the PK from the
@@ -443,6 +448,21 @@ for k=1:length(check_1),
end
end
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Checking NRDOSES and warning if used so that user can check if it was
+% done correctly.
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+% Check if NRDOSES somewhere different from NaN
+NRDOSES = data.NRDOSES(~isnan(data.NRDOSES));
+if ~isempty(NRDOSES),
+ fprintf('\nThe NRDOSES column has been defined in the dataset. Please consider expansion of these dose records to individual dose records when\n');
+ fprintf('generating the task dataset. The function IQMconvertGeneral2TaskDataset allows the definition of a flag to do this expansion automatically.\n');
+ fprintf('The expansion is beneficial since the VPC functions in IQM tools require definition of single dose records. Parameter estimation in NONMEM\n');
+ fprintf('and Monolix will not be affected by your choice. Also note that NRDOSES is equivalent to the ADDL column in NONMEM and MONOLIX and should\n');
+ fprintf('be defined accordingly (NRDOSES=N codes for 1 dose and N additional doses => N+1 total doses).\n');
+end
+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Final message
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
diff --git a/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMcheckGeneralDataFormatHeader.m b/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMcheckGeneralDataFormatHeader.m
index bef0f50cd58e3561fd809da7ad0e868a614d93ae..52cc7a572639cc3e9657a1bd9b4067ae4dafa6fb 100644
--- a/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMcheckGeneralDataFormatHeader.m
+++ b/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMcheckGeneralDataFormatHeader.m
@@ -53,8 +53,11 @@ function [data] = IQMcheckGeneralDataFormatHeader(data)
% or the value of other readouts. The values need to be in the units, defined in the UNIT column.
% Specific cases:
% - For concomitant medications the dose will be given
-% - For BLOQ values, 0 will be used
% - Severity levels for adverse events
+% - For BLOQ records: any value can be entered that is lower than the actual LLOQ.
+% It is not acceptable to set this value to �NaN� or �NA� since then no discrimination
+% can be made between �missing� and �<LLOQ�. For PK records on untransformed data �0�
+% is suggested. On log transformed data 0 should not be used but log(LLOQ/2) would be acceptable.
% Should not be populated if VALUETXT is populated
% VALUETXT string Text version of value (if available and useful)
% Character value as given in the CRF.
@@ -64,7 +67,9 @@ function [data] = IQMcheckGeneralDataFormatHeader(data)
% LLOQ numeric Lower limit of quantification of event defined by NAME
% ROUTE string Route of administration (iv, subcut, intramuscular, intraarticular, oral, inhaled, topical, rectal)
% INTERVAL numeric Interval of dosing, if single row should define multiple dosings
-% NRDOSES numeric Number of doses given with the specified interval, if single row should define multiple dosings
+% NRDOSES numeric Number of ADDITIONAL doses given within the specified interval,
+% NRDOSES=N codes for a total of N+1 doses with an interval as defined in the
+% column "INTERVAL" starting at the time defined in the dose record.
% COMMENT string Additional information for the observation/event
% For example:
% - For PK: concatenation of DMPK flag to exclude or not the PK from the
diff --git a/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMcreateGeneralDataset.m b/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMcreateGeneralDataset.m
index fc26fd915e050cbdd34c6db5c316cc74f97c2756..3becfc9c879ccd46364d2cfd5ab35facedba57c5 100644
--- a/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMcreateGeneralDataset.m
+++ b/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/IQMcreateGeneralDataset.m
@@ -1 +1 @@
-function [dataGeneral] = IQMcreateGeneralDataset(varargin)
% The IQM Tools' workflow functions assume a general dataset format that is
% independent of modeling activities and tools. This function here allows
% to generate an empty dataset with this structure. Additionally, in this
% help text here the format of this dataset is defined.
%
% [SYNTAX]
% [data] = IQMcreateGeneralDataset()
% [data] = IQMcreateGeneralDataset(NROWS)
%
% [INPUT]
% NRROWS: Number of rows in the generated (empty) dataset (default: 1)
%
% [OUTPUT]
% Empty dataset in the general dataset format
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Specification of general dataset format:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% COLUMN NAME TYPE DESCRIPTION
% -------------------------------------------------------------------------
% IXGDF numeric Column containing numeric indices for each record (1,2,3,4,5, ...)
% to allow for matching records in case of postprocessing the general dataset format
% IGNORE string Reason/comment related to exclusion of the sample/observation from the analysis
% USUBJID string Unique subject identifier
% COMPOUND string Name of the investigational compound
% STUDY string Study number
% STUDYDES string Study description
% PART string Part of study as defined per protocol (1 if only one part)
% EXTENS numeric Extension of the core study (0 if not extension, 1 if extension)
% CENTER numeric Center number
% SUBJECT string Subject identifier (within a center - typically not unique across whole dataset)
% INDNAME string Indication name
% TRTNAME string Name of actual treatment given
% TRTNAMER string Name of treatment to which individual was randomized
% VISIT numeric Visit number
% VISNAME string Visit name
% BASE numeric Flag indicating assessments at baseline
% (=0 for non-baseline, =1 for first baseline, =2 for second baseline, etc.)
% SCREEN numeric Flag indicating assessments at screening
% (=0 for non-screening, =1 for first screening, =2 for second screening, etc.)
% DATEDAY string Start date of event ('01-JUL-2015')
% DATETIME string Start time of event ('09:34')
% DURATION numeric Duration of event in same time units as TIMEUNIT
% NT numeric Planned time of event. Based on protocol, in the time unit defined in TIMEUNIT column
% TIME numeric Actual time of event relative to first administration, in the time unit defined in TIMEUNIT column
% TIMEUNIT string Unit of all numerical time definitions in the dataset ('hours','days','weeks','minutes')
% TYPENAME string Unique type of event
% NAME string Unique name for the event
% VALUE numeric Value of the event, defined by NAME. E.g., the given dose, the observed PK concentration,
% or the value of other readouts. The values need to be in the units, defined in the UNIT column.
% Specific cases:
% - For concomitant medications the dose will be given
% - For BLOQ values, 0 will be used
% - Severity levels for adverse events
% Should not be populated if VALUETXT is populated
% VALUETXT string Text version of value (if available and useful)
% Character value as given in the CRF.
% Should not be populated if VALUE is populated
% UNIT string Unit of the value reported in the VALUE column. For same event the same unit has to be used across the dataset.
% ULOQ numeric Upper limit of quantification of event defined by NAME
% LLOQ numeric Lower limit of quantification of event defined by NAME
% ROUTE string Route of administration (iv, subcut, intramuscular, intraarticular, oral, inhaled, topical, rectal)
% INTERVAL numeric Interval of dosing, if single row should define multiple dosings
% NRDOSES numeric Number of doses given with the specified interval, if single row should define multiple dosings
% COMMENT string Additional information for the observation/event
% For example:
% - For PK: concatenation of DMPK flag to exclude or not the PK from the
% DMPK analysis and comment for each sample (e.g., if the sample is
% flagged as 'Excluded' than this word should be a prefix to the comment)
% - For adverse events concatenation of the seriousness and if the AE is drug related or not
% - For an imputation, it should mention 'Imputed'.
%
% The general data format might also contain the following columns, which
% are numeric equivalents to some string columns. If not provided, then
% these can be generated automatically by using the command
% IQMconvertGeneral2TaskDataset.
%
% IND: Numeric indication flag (unique for each entry in INDNAME)
% STUDYN: Numeric study flag (unique for each entry in STUDY)
% TRT: Numeric actual treatment flag (unique for each entry in TRTNAME)
% TRTR: Numeric randomized treatment flag (unique for each entry in TRTNAMER)
% <<<COPYRIGHTSTATEMENT - IQM TOOLS PRO>>>
NROWS = 1;
if nargin==1,
NROWS = varargin{1};
end
% Create standard dataset structure
dataGeneral = table();
% Indices
dataGeneral.IXGDF = [1:NROWS]';
% General
dataGeneral.IGNORE = cell(NROWS,1); dataGeneral.IGNORE(1:end) = {''};
% Unique subject identifier
dataGeneral.USUBJID = cell(NROWS,1); dataGeneral.USUBJID(1:end) = {''};
% Study information
dataGeneral.COMPOUND = cell(NROWS,1); dataGeneral.COMPOUND(1:end) = {''};
dataGeneral.STUDY = cell(NROWS,1); dataGeneral.STUDY(1:end) = {''};
dataGeneral.STUDYDES = cell(NROWS,1); dataGeneral.STUDYDES(1:end) = {''};
dataGeneral.PART = cell(NROWS,1); dataGeneral.PART(1:end) = {''};
dataGeneral.EXTENS = zeros(NROWS,1);
dataGeneral.CENTER = NaN(NROWS,1);
dataGeneral.SUBJECT = cell(NROWS,1); dataGeneral.SUBJECT(1:end) = {''};
dataGeneral.INDNAME = cell(NROWS,1); dataGeneral.INDNAME(1:end) = {''};
% Treatment group information
dataGeneral.TRTNAME = cell(NROWS,1); dataGeneral.TRTNAME(1:end) = {''};
dataGeneral.TRTNAMER = cell(NROWS,1); dataGeneral.TRTNAMER(1:end) = {''};
% Visit information
dataGeneral.VISIT = NaN(NROWS,1);
dataGeneral.VISNAME = cell(NROWS,1); dataGeneral.VISNAME(1:end) = {''};
dataGeneral.BASE = zeros(NROWS,1);
dataGeneral.SCREEN = zeros(NROWS,1);
% Event time information
dataGeneral.DATEDAY = cell(NROWS,1); dataGeneral.DATEDAY(1:end) = {''};
dataGeneral.DATETIME = cell(NROWS,1); dataGeneral.DATETIME(1:end) = {''};
dataGeneral.DURATION = zeros(NROWS,1);
dataGeneral.NT = NaN(NROWS,1);
dataGeneral.TIME = NaN(NROWS,1);
dataGeneral.TIMEUNIT = cell(NROWS,1); dataGeneral.TIMEUNIT(1:end) = {''};
% Event value information
dataGeneral.TYPENAME = cell(NROWS,1); dataGeneral.TYPENAME(1:end) = {''};
dataGeneral.NAME = cell(NROWS,1); dataGeneral.NAME(1:end) = {''};
dataGeneral.VALUE = NaN(NROWS,1);
dataGeneral.VALUETXT = cell(NROWS,1); dataGeneral.VALUETXT(1:end) = {''};
dataGeneral.UNIT = cell(NROWS,1); dataGeneral.UNIT(1:end) = {''};
dataGeneral.ULOQ = NaN(NROWS,1);
dataGeneral.LLOQ = NaN(NROWS,1);
dataGeneral.ROUTE = cell(NROWS,1); dataGeneral.ROUTE(1:end) = {''};
dataGeneral.INTERVAL = NaN(NROWS,1);
dataGeneral.NRDOSES = NaN(NROWS,1);
% Comment
dataGeneral.COMMENT = cell(NROWS,1); dataGeneral.COMMENT(1:end) = {''};
\ No newline at end of file
+function [dataGeneral] = IQMcreateGeneralDataset(varargin)
% The IQM Tools' workflow functions assume a general dataset format that is
% independent of modeling activities and tools. This function here allows
% to generate an empty dataset with this structure. Additionally, in this
% help text here the format of this dataset is defined.
%
% [SYNTAX]
% [data] = IQMcreateGeneralDataset()
% [data] = IQMcreateGeneralDataset(NROWS)
%
% [INPUT]
% NRROWS: Number of rows in the generated (empty) dataset (default: 1)
%
% [OUTPUT]
% Empty dataset in the general dataset format
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Specification of general dataset format:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% COLUMN NAME TYPE DESCRIPTION
% -------------------------------------------------------------------------
% IXGDF numeric Column containing numeric indices for each record (1,2,3,4,5, ...)
% to allow for matching records in case of postprocessing the general dataset format
% IGNORE string Reason/comment related to exclusion of the sample/observation from the analysis
% USUBJID string Unique subject identifier
% COMPOUND string Name of the investigational compound
% STUDY string Study number
% STUDYDES string Study description
% PART string Part of study as defined per protocol (1 if only one part)
% EXTENS numeric Extension of the core study (0 if not extension, 1 if extension)
% CENTER numeric Center number
% SUBJECT string Subject identifier (within a center - typically not unique across whole dataset)
% INDNAME string Indication name
% TRTNAME string Name of actual treatment given
% TRTNAMER string Name of treatment to which individual was randomized
% VISIT numeric Visit number
% VISNAME string Visit name
% BASE numeric Flag indicating assessments at baseline
% (=0 for non-baseline, =1 for first baseline, =2 for second baseline, etc.)
% SCREEN numeric Flag indicating assessments at screening
% (=0 for non-screening, =1 for first screening, =2 for second screening, etc.)
% DATEDAY string Start date of event ('01-JUL-2015')
% DATETIME string Start time of event ('09:34')
% DURATION numeric Duration of event in same time units as TIMEUNIT
% NT numeric Planned time of event. Based on protocol, in the time unit defined in TIMEUNIT column
% TIME numeric Actual time of event relative to first administration, in the time unit defined in TIMEUNIT column
% TIMEUNIT string Unit of all numerical time definitions in the dataset ('hours','days','weeks','minutes')
% TYPENAME string Unique type of event
% NAME string Unique name for the event
% VALUE numeric Value of the event, defined by NAME. E.g., the given dose, the observed PK concentration,
% or the value of other readouts. The values need to be in the units, defined in the UNIT column.
% Specific cases:
% - For concomitant medications the dose will be given
% - Severity levels for adverse events
% - For BLOQ records: any value can be entered that is lower than the actual LLOQ.
% It is not acceptable to set this value to �NaN� or �NA� since then no discrimination
% can be made between �missing� and �<LLOQ�. For PK records on untransformed data �0�
% is suggested. On log transformed data 0 should not be used but log(LLOQ/2) would be acceptable.
% Should not be populated if VALUETXT is populated
% VALUETXT string Text version of value (if available and useful)
% Character value as given in the CRF.
% Should not be populated if VALUE is populated
% UNIT string Unit of the value reported in the VALUE column. For same event the same unit has to be used across the dataset.
% ULOQ numeric Upper limit of quantification of event defined by NAME
% LLOQ numeric Lower limit of quantification of event defined by NAME
% ROUTE string Route of administration (iv, subcut, intramuscular, intraarticular, oral, inhaled, topical, rectal)
% INTERVAL numeric Interval of dosing, if single row should define multiple dosings
% NRDOSES numeric Number of ADDITIONAL doses given within the specified interval,
% NRDOSES=N codes for a total of N+1 doses with an interval as defined in the
% column "INTERVAL" starting at the time defined in the dose record.
% COMMENT string Additional information for the observation/event
% For example:
% - For PK: concatenation of DMPK flag to exclude or not the PK from the
% DMPK analysis and comment for each sample (e.g., if the sample is
% flagged as 'Excluded' than this word should be a prefix to the comment)
% - For adverse events concatenation of the seriousness and if the AE is drug related or not
% - For an imputation, it should mention 'Imputed'.
%
% The general data format might also contain the following columns, which
% are numeric equivalents to some string columns. If not provided, then
% these can be generated automatically by using the command
% IQMconvertGeneral2TaskDataset.
%
% IND: Numeric indication flag (unique for each entry in INDNAME)
% STUDYN: Numeric study flag (unique for each entry in STUDY)
% TRT: Numeric actual treatment flag (unique for each entry in TRTNAME)
% TRTR: Numeric randomized treatment flag (unique for each entry in TRTNAMER)
% <<<COPYRIGHTSTATEMENT - IQM TOOLS PRO>>>
NROWS = 1;
if nargin==1,
NROWS = varargin{1};
end
% Create standard dataset structure
dataGeneral = table();
% Indices
dataGeneral.IXGDF = [1:NROWS]';
% General
dataGeneral.IGNORE = cell(NROWS,1); dataGeneral.IGNORE(1:end) = {''};
% Unique subject identifier
dataGeneral.USUBJID = cell(NROWS,1); dataGeneral.USUBJID(1:end) = {''};
% Study information
dataGeneral.COMPOUND = cell(NROWS,1); dataGeneral.COMPOUND(1:end) = {''};
dataGeneral.STUDY = cell(NROWS,1); dataGeneral.STUDY(1:end) = {''};
dataGeneral.STUDYDES = cell(NROWS,1); dataGeneral.STUDYDES(1:end) = {''};
dataGeneral.PART = cell(NROWS,1); dataGeneral.PART(1:end) = {''};
dataGeneral.EXTENS = zeros(NROWS,1);
dataGeneral.CENTER = NaN(NROWS,1);
dataGeneral.SUBJECT = cell(NROWS,1); dataGeneral.SUBJECT(1:end) = {''};
dataGeneral.INDNAME = cell(NROWS,1); dataGeneral.INDNAME(1:end) = {''};
% Treatment group information
dataGeneral.TRTNAME = cell(NROWS,1); dataGeneral.TRTNAME(1:end) = {''};
dataGeneral.TRTNAMER = cell(NROWS,1); dataGeneral.TRTNAMER(1:end) = {''};
% Visit information
dataGeneral.VISIT = NaN(NROWS,1);
dataGeneral.VISNAME = cell(NROWS,1); dataGeneral.VISNAME(1:end) = {''};
dataGeneral.BASE = zeros(NROWS,1);
dataGeneral.SCREEN = zeros(NROWS,1);
% Event time information
dataGeneral.DATEDAY = cell(NROWS,1); dataGeneral.DATEDAY(1:end) = {''};
dataGeneral.DATETIME = cell(NROWS,1); dataGeneral.DATETIME(1:end) = {''};
dataGeneral.DURATION = zeros(NROWS,1);
dataGeneral.NT = NaN(NROWS,1);
dataGeneral.TIME = NaN(NROWS,1);
dataGeneral.TIMEUNIT = cell(NROWS,1); dataGeneral.TIMEUNIT(1:end) = {''};
% Event value information
dataGeneral.TYPENAME = cell(NROWS,1); dataGeneral.TYPENAME(1:end) = {''};
dataGeneral.NAME = cell(NROWS,1); dataGeneral.NAME(1:end) = {''};
dataGeneral.VALUE = NaN(NROWS,1);
dataGeneral.VALUETXT = cell(NROWS,1); dataGeneral.VALUETXT(1:end) = {''};
dataGeneral.UNIT = cell(NROWS,1); dataGeneral.UNIT(1:end) = {''};
dataGeneral.ULOQ = NaN(NROWS,1);
dataGeneral.LLOQ = NaN(NROWS,1);
dataGeneral.ROUTE = cell(NROWS,1); dataGeneral.ROUTE(1:end) = {''};
dataGeneral.INTERVAL = NaN(NROWS,1);
dataGeneral.NRDOSES = NaN(NROWS,1);
% Comment
dataGeneral.COMMENT = cell(NROWS,1); dataGeneral.COMMENT(1:end) = {''};
\ No newline at end of file
diff --git a/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/auxiliary/expandGeneralNR_DOSES_intervalIQM.m b/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/auxiliary/expandGeneralNR_DOSES_intervalIQM.m
index 636a758cce4800094c03b455e6519d952cb92f22..fbabc017de44674185aebf30cb6cdc4ff85fac7f 100644
--- a/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/auxiliary/expandGeneralNR_DOSES_intervalIQM.m
+++ b/IQM Tools/IQMpro/tools/03-ClinicalData/01-GeneralDataFormat/auxiliary/expandGeneralNR_DOSES_intervalIQM.m
@@ -5,6 +5,9 @@ function [dataOut] = expandGeneralNR_DOSES_intervalIQM(data,DOSENAME)
% VIST, VISNAME, DATEDAY/TIME and DURATION entries are kept on the
% original dose definitions, as well as IXGDF. TIME and NT will be adjusted.
%
+% NOTE THAT NRDOSES codes for ADDITIONAL doses ... so NRDOSES=1 codes for 2
+% doses, etc. Spacing of doses is defined by the column INTERVAL.
+%
% [SYNTAX]
% [dataOut] = expandGeneralNR_DOSES_intervalIQM(data,DOSENAME)
%
@@ -27,16 +30,19 @@ end
dataDOSES = data(strcmp(data.NAME,DOSENAME),:);
dataOTHER = data(~strcmp(data.NAME,DOSENAME),:);
% Expand doses if both INTERVAL and NRDOSES defined
+% NRDOSES+1 DOSES need to be created if expansion is done...
dataDOSES_expanded = table();
for k=1:height(dataDOSES),
datak = dataDOSES(k,:);
if ~isnan(datak.INTERVAL) && ~isnan(datak.NRDOSES),
- dataexp = datak(ones(1,datak.NRDOSES),:);
- dataexp.TIME = dataexp.TIME+[0:datak.INTERVAL:datak.INTERVAL*(datak.NRDOSES-1)]';
- dataexp.NT = dataexp.NT+[0:datak.INTERVAL:datak.INTERVAL*(datak.NRDOSES-1)]';
+ NREXPANDED_DOSES = datak.NRDOSES + 1;
+
+ dataexp = datak(ones(1,NREXPANDED_DOSES),:);
+ dataexp.TIME = dataexp.TIME+[0:datak.INTERVAL:datak.INTERVAL*(NREXPANDED_DOSES-1)]';
+ dataexp.NT = dataexp.NT+[0:datak.INTERVAL:datak.INTERVAL*(NREXPANDED_DOSES-1)]';
dataexp.INTERVAL(1:end) = NaN;
- dataexp.NRDOSES(1:end) = NaN;
- dataDOSES_expanded = [dataDOSES_expanded; dataexp];
+ dataexp.NRDOSES(1:end) = NaN;
+ dataDOSES_expanded = [dataDOSES_expanded; dataexp];
elseif isnan(datak.INTERVAL) && isnan(datak.NRDOSES),
% Nothing to expand (single dose definition)
dataDOSES_expanded = [dataDOSES_expanded; datak];
diff --git a/IQM Tools/IQMpro/tools/03-ClinicalData/03-TaskDataset/IQMconvertGeneral2TaskDataset.m b/IQM Tools/IQMpro/tools/03-ClinicalData/03-TaskDataset/IQMconvertGeneral2TaskDataset.m
index 1dc63600882404ffe70d3d8823736fe51075432e..37e04cd8c360cc9b17fa5fa57c687a19b4f1a177 100644
--- a/IQM Tools/IQMpro/tools/03-ClinicalData/03-TaskDataset/IQMconvertGeneral2TaskDataset.m
+++ b/IQM Tools/IQMpro/tools/03-ClinicalData/03-TaskDataset/IQMconvertGeneral2TaskDataset.m
@@ -75,11 +75,11 @@ function [dataOut] = IQMconvertGeneral2TaskDataset(data,DOSENAMES,OBSNAMES,covar
% DOSE: Time dependent DOSE amount (carry forward imputation)
% This column is only present in the case that a single
% element in DOSENAMES is defined.
-% DOSE_"name" where "name" is generated from the DOSENAMES entries.
+% DOSE"name" where "name" is generated from the DOSENAMES entries.
% These columns are only present if DOSENAMES contains
% multiple entries.
%
-% For DOSE and DOSE_"name" columns. Prior to the first dose the entries are
+% For DOSE and DOSE"name" columns. Prior to the first dose the entries are
% NaN, since dose undefined.
%
% [SYNTAX]
@@ -234,7 +234,7 @@ if length(DOSENAMES)==1,
dataOut.DOSE(isnan(dataOut.DOSE)) = 0;
else
for k=1:length(DOSENAMES),
- colname_dose = sprintf('DOSE_%s',makeVariableNameIQM(DOSENAMES{k}));
+ colname_dose = sprintf('DOSE%s',makeVariableNameIQM(DOSENAMES{k}));
dataOut = IQMdataAddTimeDependentCovariate(dataOut,{DOSENAMES{k} colname_dose},FLAG_CARRY_FIRST_NON_NAN_BACKWARD);
% Handle cases where a subject has not received any dose in the dataset (DOSE values are NaN) ... set these to 0, since
diff --git a/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreBLLOQdata.m b/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreBLLOQdata.m
index 391b94b93cb4a40ada891f0a8ee6b1a084813868..a14041c19a942dcf8cb1899f65d10a8b73218912 100644
--- a/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreBLLOQdata.m
+++ b/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreBLLOQdata.m
@@ -32,7 +32,7 @@ if nargin<2,
filename = '';
end
-% Find all NAMEs which have LLOQ information
+% Find all NAMEs which have LLOQ information (only observations considered)
NAMES_BLLOQ_present = unique(data.NAME(~isnan(data.LLOQ) & data.EVID==0));
% Find all records that are BLOQ for all names that have LLOQ information
diff --git a/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreIndivData.m b/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreIndivData.m
index 3f78104f7265bee143f1574078a1b88875528324..7d95310bbca16b2a7c148294dd13886849799daf 100644
--- a/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreIndivData.m
+++ b/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreIndivData.m
@@ -193,6 +193,26 @@ if ~isempty(titlefontsize),
options.titlefontsize = titlefontsize;
end
+%% Handle case of BLOQ data present
+BLLOQdataPresent = 0;
+% Check if dataPlot contains LLOQ information and if yes then check if BLLOQ data are present
+if sum(~isnan(dataPlot.LLOQ)) > 0,
+ % LLOQ information present
+ if sum(dataPlot.DV<dataPlot.LLOQ) > 0,
+ % BLLOQ data present
+ BLLOQdataPresent = 1;
+ end
+end
+
+if BLLOQdataPresent,
+ dataPlot.isBLLOQ = double(dataPlot.DV<dataPlot.LLOQ);
+ options.nameColorGroup = 'isBLLOQ';
+ options.linecolorsCustom = [0.2 0.2 0.2; 0.8500 0.3250 0.0980];
+ options.linetypesCustom = {'.','x'};
+ options.showmarkers = 1;
+ options.markersize = 12;
+end
+
%% Do the plotting
IQMplottrellis(dataPlot,nameGroup,nameX,nameY,options)
diff --git a/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreIndivDataRelation.m b/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreIndivDataRelation.m
index 7b21f45009e98ac808d5b227a3cbc1644e32cc73..1757e08e36fce502a4b704b553a75c71acd65d46 100644
--- a/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreIndivDataRelation.m
+++ b/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreIndivDataRelation.m
@@ -3,8 +3,9 @@ function [] = IQMexploreIndivDataRelation(data,OBSNAMES,GROUP,options)
% dataset used in IQM Tools.
%
% [SYNTAX]
-% [] = IQMexploreIndivData(data,OBSNAMES,GROUP)
-% [] = IQMexploreIndivData(data,OBSNAMES,GROUP,options)
+% [] = IQMexploreIndivDataRelation(data,OBSNAMES)
+% [] = IQMexploreIndivDataRelation(data,OBSNAMES,GROUP)
+% [] = IQMexploreIndivDataRelation(data,OBSNAMES,GROUP,options)
%
% [INPUT]
% data: MATLAB PKPD dataset in general standard data format.
diff --git a/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexplorePKdataWrapper.m b/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexplorePKdataWrapper.m
index c1c9cf4ca6ce50e1f1874887d5166b0c6f69c430..201ea4ffeb2d7b00415b30a213b5d5a8fbb683fc 100644
--- a/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexplorePKdataWrapper.m
+++ b/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexplorePKdataWrapper.m
@@ -183,6 +183,7 @@ IQMexploreDataContents(data,DOSENAME,OBSNAME,[outputPath '/00_Data_Contents']);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Summary plot individual data on linear Y axis
+% Show BLLOQ data by different marker and color
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
options = [];
options.logY = 0;
@@ -197,6 +198,7 @@ IQMexploreIndivData(data,OBSNAME,'',options)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Summary plot individual data on logarithmic Y axis
+% Show BLLOQ data by different marker and color
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
options = [];
options.logY = 1;
@@ -211,6 +213,7 @@ IQMexploreIndivData(data,OBSNAME,'',options)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Plot individual data - 1 page per ID - linear Y axis
+% Show BLLOQ data by different marker and color
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
options = [];
options.logY = 0;
@@ -224,6 +227,7 @@ IQMexploreIndivData(data,OBSNAME,DOSENAME,options)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Plot individual data - 1 page per ID - log Y axis
+% Show BLLOQ data by different marker and color
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
options = [];
options.logY = 1;
@@ -236,7 +240,8 @@ options.filename = [outputPath '/04_Individual_Single_Log'];
IQMexploreIndivData(data,OBSNAME,DOSENAME,options)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-% Assessment of data availability TRT/STUDY
+%% Assessment of data availability TRT/STUDY
+% Show BLLOQ data by different marker and color
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
filename = [outputPath '/05_Summary_Study_Treatment'];
IQMstartNewPrintFigure(filename);
@@ -267,13 +272,35 @@ options.maxlegendentries = 20;
options.titlefontsize = 8;
options.labeltextsize = 8;
options.ticklabeltextsize = 8;
+
+% Handle case of BLOQ data present
+BLLOQdataPresent = 0;
+% Check if dataPlot contains LLOQ information and if yes then check if BLLOQ data are present
+if sum(~isnan(dataPlot.LLOQ)) > 0,
+ % LLOQ information present
+ if sum(dataPlot.DV<dataPlot.LLOQ) > 0,
+ % BLLOQ data present
+ BLLOQdataPresent = 1;
+ end
+end
+if BLLOQdataPresent,
+ dataPlot.isBLLOQ = double(dataPlot.DV<dataPlot.LLOQ);
+ options.nameColorGroup = 'isBLLOQ';
+ options.linecolorsCustom = [0.2 0.2 0.2; 0.8500 0.3250 0.0980];
+ options.linetypesCustom = {'.-','x-'};
+ options.showmarkers = 1;
+ options.markersize = 12;
+end
+
+% Plot
IQMplotfacetgrid(dataPlot,nameX,nameY,nameGroupX,nameGroupY,options)
IQMprintFigure(gcf,filename);
close(gcf);
IQMconvert2pdf(filename);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-% Dose normalized plots - over TIME
+%% Dose normalized plots - over TIME
+% Not considering BLOQ data (remove it prior to dose-normalization)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
filename = [outputPath '/06_Dose_Normalized_TIME'];
@@ -283,6 +310,9 @@ IQMstartNewPrintFigure(filename);
% Linear
%%%%%%%%%
dataPlot = subsetIQM(data,'NAME',OBSNAME);
+% Remove BLLOQ data
+dataPlot(dataPlot.DV<dataPlot.LLOQ,:) = [];
+
% Dose normalize the PK data
DOSE = dataPlot.DOSE; DOSE(DOSE==0) = 1;
dataPlot.DVnorm = dataPlot.DV./DOSE;
@@ -299,7 +329,7 @@ options.linewidth = 1;
options.nameSubGroup = 'USUBJID';
options.nameColorGroup = 'TRTNAME';
options.xlabelText = sprintf('Time [%s]',dataPlot.TIMEUNIT{1});
-options.ylabelText = sprintf('(DOSE NORMALIZED) %s',dataPlot.NAME{1});
+options.ylabelText = sprintf('(DOSE NORMALIZED) %s - BLLOQ data not considered',dataPlot.NAME{1});
options.ylabelfirstonly = 1;
options.logY = 0;
options.showmedian = 1;
@@ -323,6 +353,9 @@ close(gcf);
% Log
%%%%%%%%%
dataPlot = subsetIQM(data,'NAME',OBSNAME);
+% Remove BLLOQ data
+dataPlot(dataPlot.DV<dataPlot.LLOQ,:) = [];
+
% Dose normalize the PK data
DOSE = dataPlot.DOSE; DOSE(DOSE==0) = 1;
dataPlot.DVnorm = dataPlot.DV./DOSE;
@@ -339,7 +372,7 @@ options.linewidth = 1;
options.nameSubGroup = 'USUBJID';
options.nameColorGroup = 'TRTNAME';
options.xlabelText = sprintf('Time [%s]',dataPlot.TIMEUNIT{1});
-options.ylabelText = sprintf('(DOSE NORMALIZED) %s',dataPlot.NAME{1});
+options.ylabelText = sprintf('(DOSE NORMALIZED) %s - BLLOQ data not considered',dataPlot.NAME{1});
options.ylabelfirstonly = 1;
options.logY = 1;
options.showmedian = 1;
@@ -363,7 +396,8 @@ IQMconvert2pdf(filename);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-% Dose normalized plots - over TAD
+%% Dose normalized plots - over TAD
+% Not considering BLOQ data (remove it prior to dose-normalization)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
filename = [outputPath '/07_Dose_Normalized_TAD'];
@@ -373,6 +407,9 @@ IQMstartNewPrintFigure(filename);
% Linear
%%%%%%%%%
dataPlot = subsetIQM(data,'NAME',OBSNAME);
+% Remove BLLOQ data
+dataPlot(dataPlot.DV<dataPlot.LLOQ,:) = [];
+
% Dose normalize the PK data
DOSE = dataPlot.DOSE; DOSE(DOSE==0) = 1;
dataPlot.DVnorm = dataPlot.DV./DOSE;
@@ -389,7 +426,7 @@ options.linewidth = 1;
options.nameSubGroup = 'USUBJID';
options.nameColorGroup = 'TRTNAME';
options.xlabelText = sprintf('TAD [%s]',dataPlot.TIMEUNIT{1});
-options.ylabelText = sprintf('(DOSE NORMALIZED) %s',dataPlot.NAME{1});
+options.ylabelText = sprintf('(DOSE NORMALIZED) %s - BLLOQ data not considered',dataPlot.NAME{1});
options.ylabelfirstonly = 1;
options.logY = 0;
options.showmedian = 1;
@@ -413,6 +450,9 @@ close(gcf);
% Log
%%%%%%%%%
dataPlot = subsetIQM(data,'NAME',OBSNAME);
+% Remove BLLOQ data
+dataPlot(dataPlot.DV<dataPlot.LLOQ,:) = [];
+
% Dose normalize the PK data
DOSE = dataPlot.DOSE; DOSE(DOSE==0) = 1;
dataPlot.DVnorm = dataPlot.DV./DOSE;
@@ -429,7 +469,7 @@ options.linewidth = 1;
options.nameSubGroup = 'USUBJID';
options.nameColorGroup = 'TRTNAME';
options.xlabelText = sprintf('TAD [%s]',dataPlot.TIMEUNIT{1});
-options.ylabelText = sprintf('(DOSE NORMALIZED) %s',dataPlot.NAME{1});
+options.ylabelText = sprintf('(DOSE NORMALIZED) %s - BLLOQ data not considered',dataPlot.NAME{1});
options.ylabelfirstonly = 1;
options.logY = 1;
options.showmedian = 1;
@@ -452,7 +492,7 @@ close(gcf);
IQMconvert2pdf(filename);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-% Summary statistics covariates
+%% Summary statistics covariates
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if ~isempty(covNames) || ~isempty(catNames),
filename = [outputPath '/08_Summary_Statistics'];
@@ -470,9 +510,13 @@ end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Graphical exploration of potential covariate effect on PK
% We only look at dose normalized information
+% Not considering BLOQ data (remove it prior to dose-normalization)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
dataPlot = subsetIQM(data,'NAME',OBSNAME);
+% Remove BLLOQ data
+dataPlot(dataPlot.DV<dataPlot.LLOQ,:) = [];
+
% Dose normalize the PK data
DOSE = dataPlot.DOSE; DOSE(DOSE==0) = 1;
dataPlot.DVnorm = dataPlot.DV./DOSE;
@@ -520,7 +564,7 @@ if ~isempty(covNames),
options.nameSubGroup = 'USUBJID';
options.nameColorGroup = newCatCovNames{k};
options.xlabelText = sprintf('TAD [%s]',dataPlotCov.TIMEUNIT{1});
- options.ylabelText = sprintf('(DOSE NORMALIZED) %s',dataPlotCov.NAME{1});
+ options.ylabelText = sprintf('(DOSE NORMALIZED) %s - BLLOQ data not considered',dataPlotCov.NAME{1});
options.ylabelfirstonly = 1;
options.logY = 0;
options.showmedian = 1;
@@ -549,7 +593,7 @@ if ~isempty(covNames),
options.nameSubGroup = 'ID';
options.nameColorGroup = newCatCovNames{k};
options.xlabelText = sprintf('TAD [%s]',dataPlotCov.TIMEUNIT{1});
- options.ylabelText = sprintf('(DOSE NORMALIZED) %s',dataPlotCov.NAME{1});
+ options.ylabelText = sprintf('(DOSE NORMALIZED) %s - BLLOQ data not considered',dataPlotCov.NAME{1});
options.ylabelfirstonly = 1;
options.logY = 1;
options.showmedian = 1;
@@ -598,7 +642,7 @@ if ~isempty(catNames),
options.nameSubGroup = 'ID';
options.nameColorGroup = catNames{k};
options.xlabelText = sprintf('TAD [%s]',dataPlotCov.TIMEUNIT{1});
- options.ylabelText = sprintf('(DOSE NORMALIZED) %s',dataPlotCov.NAME{1});
+ options.ylabelText = sprintf('(DOSE NORMALIZED) %s - BLLOQ data not considered',dataPlotCov.NAME{1});
options.ylabelfirstonly = 1;
options.logY = 0;
options.showmedian = 1;
@@ -627,7 +671,7 @@ if ~isempty(catNames),
options.nameSubGroup = 'ID';
options.nameColorGroup = catNames{k};
options.xlabelText = sprintf('TAD [%s]',dataPlotCov.TIMEUNIT{1});
- options.ylabelText = sprintf('(DOSE NORMALIZED) %s',dataPlotCov.NAME{1});
+ options.ylabelText = sprintf('(DOSE NORMALIZED) %s - BLLOQ data not considered',dataPlotCov.NAME{1});
options.ylabelfirstonly = 1;
options.logY = 1;
options.showmedian = 1;
diff --git a/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreSummaryStats.m b/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreSummaryStats.m
index 26fec432509acf731c6692d8ffae812e0323e4a3..521a0ff40584b61459888804424529c53808a22f 100644
--- a/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreSummaryStats.m
+++ b/IQM Tools/IQMpro/tools/03-ClinicalData/04-DataExploration/IQMexploreSummaryStats.m
@@ -1,8 +1,10 @@
function [continuousTable,categoricalTable] = IQMexploreSummaryStats(data,covNames,catNames,filename)
% This function produces summary statistics for the provided dataset
% and displays the results in a table in the MATLAB window. If a filename
-% is provdided, the results are also exported to this file. The data need
-% to be provided at least in the general dataset format or in the task
+% is provdided, the results are also exported to separate files for
+% continuous and categorical covariates. (_continuous and _categorical are
+% postfixed to the filename).
+% The data need to be provided at least in the general dataset format or in the task
% specific dataset format. The covariates considered need to be available
% as columns.
%
@@ -141,11 +143,13 @@ disp(textDisplay);
textDisplay = IQMconvertCellTable2ReportTable(categoricalTable,'text');
disp(textDisplay);
-% Convert to report text and export to file if filename defined
+% Convert to report text
text1 = IQMconvertCellTable2ReportTable(continuousTable,'report');
text2 = IQMconvertCellTable2ReportTable(categoricalTable,'report');
-text = sprintf('%s\r\n\r\n%s',text1,text2);
-IQMwriteText2File(text,[strrep(filename,'.txt','') '.txt']);
+
+% Export to file if filename defined
+IQMwriteText2File(text1,[strrep(filename,'.txt','') '_continuous' '.txt']);
+IQMwriteText2File(text2,[strrep(filename,'.txt','') '_categorical' '.txt']);
diff --git a/IQM Tools/IQMpro/tools/04-NLMEdataPreparation/02-BLOQhandling/IQMhandleBLOQdata.m b/IQM Tools/IQMpro/tools/04-NLMEdataPreparation/02-BLOQhandling/IQMhandleBLOQdata.m
index 6d28042618d211d536805377967b380cdbbb45bd..7e4b29e418eacd28ff4a8b0dea822c5342d6c699 100644
--- a/IQM Tools/IQMpro/tools/04-NLMEdataPreparation/02-BLOQhandling/IQMhandleBLOQdata.m
+++ b/IQM Tools/IQMpro/tools/04-NLMEdataPreparation/02-BLOQhandling/IQMhandleBLOQdata.m
@@ -63,7 +63,7 @@ end
% Copy dataset
datanew = data;
-% Find all NAMEs which have LLOQ information
+% Find all NAMEs which have LLOQ information and are observations
NAMES_BLLOQ_present = unique(datanew.NAME(~isnan(datanew.LLOQ) & datanew.EVID==0));
% Find all records that are BLOQ for all names that have LLOQ information
diff --git a/IQM Tools/IQMpro/tools/04-NLMEdataPreparation/03-DataNLMEConversion/IQMconvertTask2NLMEdataset.m b/IQM Tools/IQMpro/tools/04-NLMEdataPreparation/03-DataNLMEConversion/IQMconvertTask2NLMEdataset.m
index 875aef45f720be05a3a1220d636223947cc8c08f..c1c79039d0126fc05c9ba6072e28670c47addec7 100644
--- a/IQM Tools/IQMpro/tools/04-NLMEdataPreparation/03-DataNLMEConversion/IQMconvertTask2NLMEdataset.m
+++ b/IQM Tools/IQMpro/tools/04-NLMEdataPreparation/03-DataNLMEConversion/IQMconvertTask2NLMEdataset.m
@@ -29,7 +29,8 @@ function [dataNLME] = IQMconvertTask2NLMEdataset(data,DOSENAMES,OBSNAMES,covName
% with a more informed dataset
% - Displaying mapping between ADM and ROUTE, and YTYPE and OBSNAMES
% - If NRDOSES and INTERVAL not all NaN then an ADDL and II column are
-% added to the dataNLME dataset
+% added to the dataNLME dataset. NRDOSES is same as ADDL and INTERVAL
+% same as II.
%
% [SYNTAX]
% [dataNLME] = IQMconvertTask2NLMEdataset(data,DOSENAMES,OBSNAMES,covNames,catNames,regressionNames)
@@ -130,7 +131,7 @@ if sum(abs(dataNLME.CENS)) ~= 0,
end
% Define initial structure of NLME dataset
-varNames = {'IXGDF' 'IGNORE' 'ID' 'USUBJID' 'STUDY' 'STUDYN' 'TRT' 'TRTNAME' 'TIME' 'TIMEPOS' 'NT' 'TIMEUNIT' 'YTYPE' 'NAME' 'DV' 'UNIT' 'CENS' 'MDV' 'EVID' 'AMT' 'ADM' 'INTERVAL' 'NRDOSES' 'ROUTE' 'RATE' 'TINF'};
+varNames = {'IXGDF' 'IGNORE' 'ID' 'USUBJID' 'STUDY' 'STUDYN' 'TRT' 'TRTNAME' 'TIME' 'TIMEPOS' 'NT' 'TAD' 'TIMEUNIT' 'YTYPE' 'NAME' 'DV' 'UNIT' 'CENS' 'MDV' 'EVID' 'AMT' 'ADM' 'INTERVAL' 'NRDOSES' 'ROUTE' 'RATE' 'TINF'};
% Check if DOSE is present in the dataset then add it
if ~isempty(strmatchIQM('DOSE',dataNLME.Properties.VariableNames,'exact')),
@@ -139,7 +140,7 @@ if ~isempty(strmatchIQM('DOSE',dataNLME.Properties.VariableNames,'exact')),
else
% DOSE not present ... might be due to multiple DOSENAMES ... check
for k=1:length(DOSENAMES),
- doseNameTest = ['DOSE_' makeVariableNameIQM(DOSENAMES{k})];
+ doseNameTest = ['DOSE' makeVariableNameIQM(DOSENAMES{k})];
if ~isempty(strmatchIQM(doseNameTest,dataNLME.Properties.VariableNames,'exact')),
% if only one DOSENAME(S) then rename to DOSE ... otherwise
% keep changed name
diff --git a/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/01-MONOLIX/IQMcreateMLXTRANfile.m b/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/01-MONOLIX/IQMcreateMLXTRANfile.m
index 8c7ba80d95750ba99fde522863616d764b905705..2f4ce57625823c8d6c2c7933eef2cc82e08e5b5e 100644
--- a/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/01-MONOLIX/IQMcreateMLXTRANfile.m
+++ b/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/01-MONOLIX/IQMcreateMLXTRANfile.m
@@ -9,7 +9,7 @@ function [filename, moddosinfo] = IQMcreateMLXTRANfile(model,dosing,varargin)
% [SYNTAX]
% [filename] = IQMcreateMLXTRANfile(model,dosing)
% [filename] = IQMcreateMLXTRANfile(model,dosing,filename)
-% [filename] = IQMcreateMLXTRANfile(model,dosing,filename,SILENT,regressionParameters)
+% [filename] = IQMcreateMLXTRANfile(model,dosing,filename,SILENT,regressionParameters,startTime)
%
% [INPUT]
% model: IQMmodel (annotated with additional information, see above)
@@ -29,6 +29,8 @@ function [filename, moddosinfo] = IQMcreateMLXTRANfile(model,dosing,varargin)
% in the dataset. Default: {}. If this input argument
% is given, then the regressor information in the
% IQMmodel and IQMdosing scheme is ignored!!!
+% startTime: Allows to set the start of the integration to a desired time.
+% Default: []: do not set.
%
% [OUTPUT]
% filename: constructed from the model name "modelname_MLXTRAN".
@@ -81,17 +83,22 @@ if nargin >= 3,
end
SILENT = 0;
-if nargin == 4,
+if nargin >= 4,
SILENT = varargin{2};
end
regressionParameters = {};
IGNORE_MODELDEFINED_REGRESSORS = 0;
-if nargin == 5,
+if nargin >= 5,
regressionParameters = varargin{3};
IGNORE_MODELDEFINED_REGRESSORS = 1;
end
+startTime = [];
+if nargin>=6,
+ startTime = varargin{4};
+end
+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Start conversion protocol with important information
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@@ -392,6 +399,12 @@ fprintf(fid,'; =============================================\r\n');
fprintf(fid,'EQUATION:\r\n');
fprintf(fid,'; =============================================\r\n');
+if ~isempty(startTime),
+ fprintf(fid,'\r\n\t; Start time of integration');
+ fprintf(fid,'\r\n\t; -------------------------\r\n');
+ fprintf(fid,'\tt0 = %g\r\n\r\n',startTime);
+end
+
fprintf(fid,'\r\n\t; Always use stiff solver');
fprintf(fid,'\r\n\t; -----------------------\r\n');
fprintf(fid,'\todeType = stiff\r\n\r\n');
diff --git a/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/01-MONOLIX/IQMcreateMONOLIXproject.m b/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/01-MONOLIX/IQMcreateMONOLIXproject.m
index 81270ff061a4069495bb84956fe2e049b059ddbd..ce128ba061d4e4e8826192eba11bbbbd48f14bc1 100644
--- a/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/01-MONOLIX/IQMcreateMONOLIXproject.m
+++ b/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/01-MONOLIX/IQMcreateMONOLIXproject.m
@@ -81,6 +81,7 @@ function IQMcreateMONOLIXproject(model,dosing,data,projectPath,varargin)
% options.algorithm.LLsetting: 'linearization' (default) or 'importantsampling'
% options.algorithm.FIMsetting: 'linearization' (default) or 'stochasticApproximation'
% options.algorithm.INDIVparametersetting: 'conditionalMode' (default) ... others not considered for now.
+% options.algorithm.startTime: start time of integration. default: [] (not set).
% options.SILENT: =0: do some output in the command window, =1: do no output in command window (default: 0)
% options.keepProjectFolder: =0: remover already existing folder, =1: keep it
%
@@ -159,6 +160,7 @@ try SILENT = options.SILENT;
try INDIVparametersetting = options.algorithm.INDIVparametersetting; catch, INDIVparametersetting = 'conditionalMode'; end
try LLsetting = options.algorithm.LLsetting; catch, LLsetting = 'linearization'; end
try FIMsetting = options.algorithm.FIMsetting; catch, FIMsetting = 'linearization'; end
+try startTime = options.algorithm.startTime; catch, startTime = []; end
try keepProjectFolder = options.keepProjectFolder; catch, keepProjectFolder = 0; end
try Ntests = options.Ntests; catch, Ntests = 1; end
try std_noise_setting = options.std_noise_setting; catch, std_noise_setting = 0.5; options.std_noise_setting = 0.5; end
@@ -372,7 +374,7 @@ ms = struct(model);
% Provide the correct regression parameter names as ordered in the dataset
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
regressionParametersData = dataCSV.Properties.VariableNames(strmatchIQM('X',explodePCIQM(data.dataHeaderIdent),'exact'));
-[modelFileName, modelInfo] = IQMcreateMLXTRANfile(model,dosing,'model',SILENT,regressionParametersData);
+[modelFileName, modelInfo] = IQMcreateMLXTRANfile(model,dosing,'model',SILENT,regressionParametersData,startTime);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% If needed, reorder parameters
diff --git a/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/01-MONOLIX/auxiliary/sampleMONOLIXpopulationParametersIQM.m b/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/01-MONOLIX/auxiliary/sampleMONOLIXpopulationParametersIQM.m
index 576628fdb6faf9357f5a1e4bf17d1293e11d87b8..8cdc03c5a9e637362687ba534ce4d8077aa05025 100644
--- a/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/01-MONOLIX/auxiliary/sampleMONOLIXpopulationParametersIQM.m
+++ b/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/01-MONOLIX/auxiliary/sampleMONOLIXpopulationParametersIQM.m
@@ -343,13 +343,20 @@ COV2informationCategories = {};
COV2informationValues = {};
for k=1:length(COVparameters),
- if k==1 || isempty(strmatchIQM(COVparameters{k},COV2parameters,'exact')),
+ if k==1,
COV2parameters{k} = COVparameters{k};
COV2covNames{k} = COVcovNames{k};
% Determine reference value
reference = input.covariates.catReference(strmatchIQM(COVcovNames{k},input.covariates.catNames,'exact'));
COV2informationCategories{k} = [reference COVcategory(k)];
COV2informationValues{k} = [0 COVvalues(k)];
+ elseif isempty(strmatchIQM(COVparameters{k},COV2parameters,'exact')),
+ COV2parameters{end+1} = COVparameters{k};
+ COV2covNames{end+1} = COVcovNames{k};
+ % Determine reference value
+ reference = input.covariates.catReference(strmatchIQM(COVcovNames{k},input.covariates.catNames,'exact'));
+ COV2informationCategories{end+1} = [reference COVcategory(k)];
+ COV2informationValues{end+1} = [0 COVvalues(k)];
else
% No all parameter names are present ... need to check if covariate
% name already present in this parameter
diff --git a/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/02-NONMEM/auxiliary/modelconversion/checkAndChangeModelSyntax4NONMEMconversionIQM.m b/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/02-NONMEM/auxiliary/modelconversion/checkAndChangeModelSyntax4NONMEMconversionIQM.m
index 067acea1e397ed454608240cbef00665c15185a5..6b8008f47386b5d983a11801aa8adee255f799d7 100644
--- a/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/02-NONMEM/auxiliary/modelconversion/checkAndChangeModelSyntax4NONMEMconversionIQM.m
+++ b/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/02-NONMEM/auxiliary/modelconversion/checkAndChangeModelSyntax4NONMEMconversionIQM.m
@@ -103,9 +103,12 @@ end
available_state_numbers = setdiff(1:length(ms.states),input_numbers);
% Distribute the available_state_numbers to the NaN input values
input_states(isnan(input_states(:,1)),1) = available_state_numbers;
+
+
% Column one in input_states defines the new state numbers and second column the old ones
% Now switch the ODEs
-ms.states(input_states(:,2)) = ms.states(input_states(:,1));
+sort_ix = sortrows(input_states,1);
+ms.states = ms.states(sort_ix(:,2));
% Need to reassign stateindex in ms.states.input
for k=1:length(ms.inputs),
ms.inputs(k).stateindex = input_numbers(k);
diff --git a/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/02-NONMEM/auxiliary/parseNONMEMresultsIQM.m b/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/02-NONMEM/auxiliary/parseNONMEMresultsIQM.m
index c5b2e143b6735acac30370d02436da0c7e567e5f..f5a38fcef165d371b3c9505f2d5c768a84156806 100644
--- a/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/02-NONMEM/auxiliary/parseNONMEMresultsIQM.m
+++ b/IQM Tools/IQMpro/tools/05-NLMEtoolsInterfaces/02-NONMEM/auxiliary/parseNONMEMresultsIQM.m
@@ -28,6 +28,7 @@ function [ output ] = parseNONMEMresultsIQM( path_to_nonmem_project_folder, tran
% <<<COPYRIGHTSTATEMENT - IQM TOOLS PRO>>>
+% Set seed to lead to consistent sampling results
if nargin == 1,
transformFlag = 0;
end
diff --git a/IQM Tools/IQMpro/tools/06-NMLEprojectHandling/IQMcreateNLMEproject.m b/IQM Tools/IQMpro/tools/06-NMLEprojectHandling/IQMcreateNLMEproject.m
index f289856335b4f7bb998d79dcee5db418d2683719..ccfaeabf5ce90db98ae5b957dc0fe5d84f6ad1b5 100644
--- a/IQM Tools/IQMpro/tools/06-NMLEprojectHandling/IQMcreateNLMEproject.m
+++ b/IQM Tools/IQMpro/tools/06-NMLEprojectHandling/IQMcreateNLMEproject.m
@@ -101,6 +101,7 @@ function [] = IQMcreateNLMEproject(TOOL,model,dosing,data,projectPath,options)
% options.algorithm.LLsetting: 'linearization' (default) or 'importantsampling'
% options.algorithm.FIMsetting: 'linearization' (default) or 'stochasticApproximation'
% options.algorithm.INDIVparametersetting: 'conditionalMode' (default) ... others not considered for now.
+% options.algorithm.startTime: start time of integration. default: [] (not set).
% <<<COPYRIGHTSTATEMENT - IQM TOOLS PRO>>>
diff --git a/IQM Tools/IQMpro/tools/06-NMLEprojectHandling/IQMgetNLMEfitIndivPopMeanParam.m b/IQM Tools/IQMpro/tools/06-NMLEprojectHandling/IQMgetNLMEfitIndivPopMeanParam.m
index 52f15a2d875bfc5671e5acea4de3f03c87101381..5d6396d7a6acee5710c4bd02cc5e0031d560657a 100644
--- a/IQM Tools/IQMpro/tools/06-NMLEprojectHandling/IQMgetNLMEfitIndivPopMeanParam.m
+++ b/IQM Tools/IQMpro/tools/06-NMLEprojectHandling/IQMgetNLMEfitIndivPopMeanParam.m
@@ -47,31 +47,61 @@ elseif isempty(data),
cd(oldpath);
end
-% Get data with only covariates covNames and catNames and USUBJID
-dataCOV = unique(data(:,{'USUBJID',header.COVARIATENAMES{:}}));
+% process some things
+x = strmatchIQM('',header.CATNAMES);
+header.CATNAMES(x) = [];
+x = strmatchIQM('',header.COVNAMES);
+header.COVNAMES(x) = [];
-% Cycle through dataCOV and keep only first entry for a single USUBJID
+% Cycle through data and keep only first entry for a single USUBJID
% (should not happen anyway)
-allID = unique(dataCOV.USUBJID);
+dataFIRST = data;
+allID = unique(dataFIRST.USUBJID);
for k=1:length(allID),
- ix = ixdataIQM(dataCOV,'USUBJID',allID(k));
+ ix = ixdataIQM(dataFIRST,'USUBJID',allID(k));
if length(ix)>1,
- dataCOV(ix(2:end),:) = [];
+ dataFIRST(ix(2:end),:) = [];
end
end
-% Separate dataCOV into continuous and categorical covariates
-dataCAT = dataCOV(:,{'USUBJID' header.CATNAMES{:}});
-dataCOV = dataCOV(:,{'USUBJID' header.COVNAMES{:}});
+% Separate dataFIRST into continuous and categorical covariates
+if ~isempty(header.CATNAMES),
+ if ~isempty(header.CATNAMES{1}),
+ dataCAT = dataFIRST(:,{'USUBJID' header.CATNAMES{:}});
+ end
+else
+ dataCAT = [];
+end
+if ~isempty(header.COVNAMES),
+ if ~isempty(header.COVNAMES{1}),
+ dataCOV = dataFIRST(:,{'USUBJID' header.COVNAMES{:}});
+ end
+else
+ dataCOV = [];
+end
% Sample population parameters from fit
FLAG_SAMPLE = 3; % Use point estimates of population parameters (do not consider uncertainty)
% Return Nsamples sets of population parameters with covariates taken into account.
-NSAMPLES = height(dataCOV);
-covNames = dataCOV.Properties.VariableNames(2:end);
-covValues = table2array(dataCOV(:,2:end));
-catNames = dataCAT.Properties.VariableNames(2:end);
-catValues = table2array(dataCAT(:,2:end));
+NSAMPLES = height(dataFIRST);
+
+if isempty(dataCOV),
+ covNames = {};
+ covValues = [];
+else
+ covNames = dataCOV.Properties.VariableNames(2:end);
+ covValues = table2array(dataCOV(:,2:end));
+end
+
+if isempty(dataCAT),
+ catNames = {};
+ catValues = [];
+else
+ catNames = dataCAT.Properties.VariableNames(2:end);
+ catValues = table2array(dataCAT(:,2:end));
+end
+
+% Sample parameters
param = IQMsampleNLMEfitParam(projectPath,FLAG_SAMPLE,NSAMPLES,covNames,covValues,catNames,catValues);
% Create table out of parameters and names
@@ -79,7 +109,7 @@ popmean_param = array2table(param.parameterValuesPopulation);
popmean_param.Properties.VariableNames = param.parameterNames;
% Add USUBJID to table (same order as for covariates)
-popmean_param.USUBJID = dataCOV.USUBJID;
+popmean_param.USUBJID = dataFIRST.USUBJID;
% Reorder to have USUBJID in the front
popmean_param = popmean_param(:,[end 1:end-1]);
diff --git a/IQM Tools/IQMpro/tools/07-NLMEfitAnalysis/IQMfitanalysisConditionNumber.m b/IQM Tools/IQMpro/tools/07-NLMEfitAnalysis/IQMfitanalysisConditionNumber.m
new file mode 100644
index 0000000000000000000000000000000000000000..7e726fca419220b1b0c5b8a66b680a7ac1ccb524
--- /dev/null
+++ b/IQM Tools/IQMpro/tools/07-NLMEfitAnalysis/IQMfitanalysisConditionNumber.m
@@ -0,0 +1,78 @@
+function [cNfull,cNfebeta,cNreerror,cNcorr] = IQMfitanalysisConditionNumber(projectPath)
+% Determine the condition number for an NLME fit.
+% Several condition numbers will be generated:
+% - condition number for the full covariance matrix (cNfull)
+% - condition number for fixed effects + covariate coefficients (cNfebeta)
+% - condition number for random effects + error model parameters (cNreerror)
+% - condition number for correlation parameters (cNcorr)
+%
+% [SYNTAX]
+% [cNfull,cNfebeta,cNreerror,cNcorr] = IQMfitanalysisConditionNumber(projectPath)
+%
+% [INPUT]
+% projectPath: Path to the project.nmctl NONMEM project file
+%
+% [OUTPUT]
+% Different condition numbers
+
+% <<<COPYRIGHTSTATEMENT - IQM TOOLS PRO>>>
+
+% parse the results
+if isNONMEMprojectIQM(projectPath),
+ x = parseNONMEMresultsIQM(projectPath);
+elseif isMONOLIXprojectIQM(projectPath),
+ x = parseMONOLIXresultsIQM(projectPath);
+else
+ error('Unknown project type.');
+end
+
+cNfull = NaN;
+cNfebeta = NaN;
+cNreerror = NaN;
+cNcorr = NaN;
+
+if ~isempty(x.parameters.correlationmatrix),
+
+ cNfull = cond(x.parameters.correlationmatrix);
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % fE and beta
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ names = [x.rawParameterInfo.fixedEffects.names x.rawParameterInfo.covariate.names];
+ ix_all = [];
+ for k=1:length(names),
+ ix_all(k) = strmatchIQM(names{k},x.parameters.names,'exact');
+ end
+ cNfebeta = cond(x.parameters.correlationmatrix(ix_all,ix_all));
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % omegas and error parameters
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ if isNONMEMprojectIQM(projectPath),
+ names = [strrep(x.rawParameterInfo.randomEffects.names,'omega','omega2') x.rawParameterInfo.errorParameter.names];
+ else
+ names = [x.rawParameterInfo.randomEffects.names x.rawParameterInfo.errorParameter.names];
+ end
+ ix_all = [];
+ for k=1:length(names),
+ ix_all(k) = strmatchIQM(names{k},x.parameters.names,'exact');
+ end
+ cNreerror = cond(x.parameters.correlationmatrix(ix_all,ix_all));
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % Correlation of correlations
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ if isNONMEMprojectIQM(projectPath),
+ names = [strrep(x.rawParameterInfo.correlation.names,'corr(','omega2(')];
+ else
+ names = x.rawParameterInfo.correlation.names;
+ end
+ if length(names) >= 1,
+ ix_all = [];
+ for k=1:length(names),
+ ix_all(k) = strmatchIQM(names{k},x.parameters.names);
+ end
+ cNcorr = cond(x.parameters.correlationmatrix(ix_all,ix_all));
+ end
+end
+
diff --git a/IQM Tools/IQMpro/tools/07-NLMEfitAnalysis/IQMfitanalysisShrinkage.m b/IQM Tools/IQMpro/tools/07-NLMEfitAnalysis/IQMfitanalysisShrinkage.m
new file mode 100644
index 0000000000000000000000000000000000000000..0833667f450409d2295e326efe11508e51c83259
--- /dev/null
+++ b/IQM Tools/IQMpro/tools/07-NLMEfitAnalysis/IQMfitanalysisShrinkage.m
@@ -0,0 +1,31 @@
+function [OMEGAnames,eta_shrinkage_percent] = IQMfitanalysisShrinkage(projectPath)
+% This function determines the ETA shrinkage for the random effects
+% ETA Shrinkage is defined as:
+%
+% eta_shrinkage_percent = 100*(1-std(eta))/OMEGA;
+%
+% [SYNTAX]
+% [OMEGAnames,eta_shrinkage_percent] = IQMfitanalysisShrinkage(projectPath)
+%
+% [INPUT]
+% projectPath: Path to a NONMEM or MONOLIX project folder.
+%
+% [OUTPUT]
+% OMEGAnames: Cell-array with names of all fixed effect parameters
+% eta_shrinkage_percent: Vector with shrinkage values. NaN if parameter
+% was not having a random effect.
+
+% <<<COPYRIGHTSTATEMENT - IQM TOOLS PRO>>>
+
+% Handle NONMEM/MONOLIX
+if isMONOLIXprojectIQM(projectPath),
+ [ dataeta, OMEGA, OMEGAnames ] = parseMONOLIXetasIQM( projectPath );
+elseif isNONMEMprojectIQM(projectPath),
+ [ dataeta, OMEGA, OMEGAnames ] = parseNONMEMetasIQM( projectPath );
+else
+ error('Unknown project type.');
+end
+
+% Determine shrinkage in percent
+eta_shrinkage_percent = 100*(1-std(table2array(dataeta))./OMEGA);
+
diff --git a/IQM Tools/IQMpro/tools/08-NLMEfitSummary/IQMfitsummaryTable.m b/IQM Tools/IQMpro/tools/08-NLMEfitSummary/IQMfitsummaryTable.m
index 1dbb12a119f37d8780ad63f126c70d032081a532..175c7d0bb01011d3b571b1fbbae42a8f06fa38e4 100644
--- a/IQM Tools/IQMpro/tools/08-NLMEfitSummary/IQMfitsummaryTable.m
+++ b/IQM Tools/IQMpro/tools/08-NLMEfitSummary/IQMfitsummaryTable.m
@@ -1,11 +1,11 @@
-function [] = IQMfitsummaryTable(projectPaths,filename)
+function [tableCell] = IQMfitsummaryTable(projectPaths,filename)
% This function creates a typical report-type NLME model parameter table.
% "projectPaths" defines the NLME projects that are to be reported in this
% table.
%
% [SYNTAX]
-% [] = IQMfitsummaryTable(projectPaths)
-% [] = IQMfitsummaryTable(projectPaths,filename)
+% [tableCell] = IQMfitsummaryTable(projectPaths)
+% [tableCell] = IQMfitsummaryTable(projectPaths,filename)
%
% [INPUT]
% projectPaths: Cell-array with paths to NLME projects to create the
@@ -110,7 +110,7 @@ end
tableCell = {'<TT>' 'Population parameter estimates for considered models'};
tableCell(end+1,1:2+length(RESULTS)) = {'<TH>' '' modelsNames_Short{:}};
-tableCell(end+1,1:2) = {'<TH>' 'Parameter'};
+tableCell(end+1,1:2) = {'<TR>' 'Parameter'};
for k=1:length(RESULTS),
tableCell(end,2+k) = {'Value (RSE%) (trans)'};
end
diff --git a/IQM Tools/IQMpro/tools/08-NLMEfitSummary/IQMfitsummaryTableSingle.m b/IQM Tools/IQMpro/tools/08-NLMEfitSummary/IQMfitsummaryTableSingle.m
new file mode 100644
index 0000000000000000000000000000000000000000..9bd8575bceecfac3c423f21376d8ecd36f1c633a
--- /dev/null
+++ b/IQM Tools/IQMpro/tools/08-NLMEfitSummary/IQMfitsummaryTableSingle.m
@@ -0,0 +1,420 @@
+function [tableCell] = IQMfitsummaryTableSingle(projectPath,filename)
+% This function creates a typical report-type NLME model parameter table.
+% "projectPaths" defines the NLME projects that are to be reported in this
+% table. In contrast to IQMfitsummaryTable this function creates a single
+% table with dedicated columns for value, RSE and trans.
+%
+% [SYNTAX]
+% [tableCell] = IQMfitsummaryTableSingle(projectPath)
+% [tableCell] = IQMfitsummaryTableSingle(projectPath,filename)
+%
+% [INPUT]
+% projectPath: String with path to NLME project to create the
+% fit results table for.
+% filename: Path and filename where to generate the output file.
+% default: '' (no file generated).
+%
+% [OUTPUT]
+% Table shown in command window and if desired also saved to file for
+% reporting purpose.
+
+% <<<COPYRIGHTSTATEMENT - IQM TOOLS PRO>>>
+
+% Handle variable input arguments
+if nargin<2,
+ filename = '';
+end
+
+% Handle cell-array thing
+if ischar(projectPath),
+ projectPath = {projectPath};
+end
+
+% Check projects if NLME
+for k=1:length(projectPath),
+ if ~isNLMEprojectIQM(projectPath{k}),
+ error('Path "%s" does not contain an NLME project.',projectPath{k});
+ end
+end
+
+% Load project results - do not reorder
+order = '';
+RESULTS = parseSelectedProjectFolderResultsIQM(projectPath,order);
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Process the RESULT information somewhat
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+% Determine all available fixed effect parameters in the models
+ALLfixEffectNames = {};
+for k=1:length(RESULTS),
+ if ~isempty(RESULTS(k).rawParameterInfo),
+ fek = RESULTS(k).rawParameterInfo.fixedEffects.names;
+ ALLfixEffectNames = [ALLfixEffectNames fek];
+ end
+end
+ALLfixEffectNames = unique(ALLfixEffectNames);
+
+% Determine all available random effect parameters in the models
+ALLrandomEffectNames = {};
+for k=1:length(RESULTS),
+ if ~isempty(RESULTS(k).rawParameterInfo),
+ fek = RESULTS(k).rawParameterInfo.randomEffects.names;
+ ALLrandomEffectNames = [ALLrandomEffectNames fek];
+ end
+end
+ALLrandomEffectNames = unique(ALLrandomEffectNames);
+
+% Determine all available correlation parameters in the models
+ALLcorrelationNames = {};
+for k=1:length(RESULTS),
+ if ~isempty(RESULTS(k).rawParameterInfo),
+ fek = RESULTS(k).rawParameterInfo.correlation.names;
+ ALLcorrelationNames = [ALLcorrelationNames fek];
+ end
+end
+ALLcorrelationNames = unique(ALLcorrelationNames);
+
+% Determine all available covariate parameters in the models
+ALLcovariateNames = {};
+for k=1:length(RESULTS),
+ if ~isempty(RESULTS(k).rawParameterInfo),
+ fek = RESULTS(k).rawParameterInfo.covariate.names;
+ ALLcovariateNames = [ALLcovariateNames fek];
+ end
+end
+ALLcovariateNames = unique(ALLcovariateNames);
+
+% Determine all available covariate parameters in the models
+ALLerrorNames = {};
+for k=1:length(RESULTS),
+ if ~isempty(RESULTS(k).rawParameterInfo),
+ fek = RESULTS(k).rawParameterInfo.errorParameter.names;
+ ALLerrorNames = [ALLerrorNames fek];
+ end
+end
+ALLerrorNames = unique(ALLerrorNames);
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Calculate shrinkage
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+[shOMEGAnames,sheta_shrinkage_percent] = IQMfitanalysisShrinkage(projectPath{1});
+
+% Need to order sheta_shrinkage_percent according to ALLrandomEffectNames
+eta_shrinkage_ordered = [];
+for k=1:length(ALLrandomEffectNames),
+ xxx = strrep(strrep(ALLrandomEffectNames{k},')',''),'omega(','');
+ % find index in ALLrandomEffectNames
+ ix = strmatchIQM(xxx,shOMEGAnames,'exact');
+ if isempty(ix),
+ error('Check!');
+ end
+ eta_shrinkage_ordered(k) = sheta_shrinkage_percent(ix);
+end
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Create the table
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Retain last name of model only
+modelsNames_Short = {RESULTS.model};
+for k=1:length(modelsNames_Short),
+ [~,f,e] = fileparts(modelsNames_Short{k});
+ modelsNames_Short{k} = [f e];
+end
+
+
+
+tableCell = {'<TT>' 'Population parameter estimates for considered model' '' '' '' ''};
+tableCell(end+1,:) = {'<TH>' 'Parameter' 'Value' 'RSE (%)' 'Transformation' 'Shrinkage (%)'};
+
+
+% FIXED EFFECTs - non-estimated get postfix (FIX).
+FLAG_FOOTER = 0;
+for k0=1:length(ALLfixEffectNames),
+ tableCell(end+1,1:2) = {'<TR>' ALLfixEffectNames{k0}};
+ for k=1:length(RESULTS),
+ % Check if parameter in model
+ ix = strmatchIQM(ALLfixEffectNames{k0},RESULTS(k).rawParameterInfo.fixedEffects.names,'exact');
+ if ~isempty(ix),
+ modelVALUES = RESULTS(k).rawParameterInfo.fixedEffects.values(ix);
+ modelRSES = RESULTS(k).rawParameterInfo.fixedEffects.rse(ix);
+ modelESTIMATED = RESULTS(k).rawParameterInfo.fixedEffects.estimated(ix);
+ if modelESTIMATED,
+ value = sprintf('%1.3g',modelVALUES);
+ if RESULTS(k).NONMEM,
+ rse = sprintf('%1.3g%%*',modelRSES);
+ FLAG_FOOTER = 1;
+ else
+ rse = sprintf('%1.3g%%',modelRSES);
+ end
+ else
+ if ~isnan(modelVALUES),
+ value = sprintf('%1.3g',modelVALUES);
+ rse = sprintf('(FIX)');
+ valuerse = sprintf('%s %s',postFillCharIQM(value,10,' '),postFillCharIQM(rse,14,' '));
+ else
+ value = NaN;
+ rse = NaN;
+ end
+ end
+ % Get MU referencing transformation of fixed effect
+ modelTRANS = RESULTS(k).rawParameterInfo.fixedEffects.distribution_info{ix};
+ if strcmp(modelTRANS,'(psi)'),
+ modelTRANS = 'normal';
+ elseif strcmp(modelTRANS,'log(psi./(1-psi))'),
+ modelTRANS = 'logitnormal';
+ elseif strcmp(modelTRANS,'log(psi)'),
+ modelTRANS = 'lognormal';
+ else
+ error('Unknown transformation of fixed effect.');
+ end
+ else
+ valuerse = '-';
+ value = '-';
+ rse = '-';
+ modelTRANS = '-';
+ end
+ % Add to table
+ tableCell(end,2+k) = {value};
+ tableCell(end,2+k+1) = {rse};
+ tableCell(end,2+k+2) = {modelTRANS};
+ end
+end
+
+% % Add footer if needed
+% if FLAG_FOOTER,
+% tableCell(end+1,1:2) = {'<TF>' '* Standard errors for NONMEM models approximated by sampling due to MU-Referencing'};
+% end
+
+% Separator
+tableCell(end+1,1) = {'<TR>'};
+
+% RANDOM EFFECTs - non-estimated get postfix (FIX).
+for k0=1:length(ALLrandomEffectNames),
+ tableCell(end+1,1:2) = {'<TR>' ALLrandomEffectNames{k0}};
+ for k=1:length(RESULTS),
+ % Check if parameter in model
+ ix = strmatchIQM(ALLrandomEffectNames{k0},RESULTS(k).rawParameterInfo.randomEffects.names,'exact');
+ if ~isempty(ix),
+ modelVALUES = RESULTS(k).rawParameterInfo.randomEffects.values(ix);
+ modelRSES = RESULTS(k).rawParameterInfo.randomEffects.rse(ix);
+ modelESTIMATED = RESULTS(k).rawParameterInfo.randomEffects.estimated(ix);
+ if modelESTIMATED==1,
+ value = sprintf('%1.3g',modelVALUES);
+ rse = sprintf('%1.3g%%',modelRSES);
+ valuerse = sprintf('%s %s',postFillCharIQM(value,10,' '),rse);
+ else
+ if ~isnan(modelVALUES),
+ value = sprintf('%1.3g',modelVALUES);
+ rse = sprintf('(FIX)',modelRSES);
+ valuerse = sprintf('%s %s',postFillCharIQM(value,10,' '),rse);
+ else
+ value = NaN;
+ rse = NaN;
+ end
+ end
+ else
+ value = '-';
+ rse = '-';
+ end
+ % Add to table
+ tableCell(end,2+k) = {value};
+ tableCell(end,2+k+1) = {rse};
+ % Add shrinkage
+ if isnan(eta_shrinkage_ordered(k0)),
+ tableCell(end,2+k+3) = {'-'};
+ else
+ tableCell{end,2+k+3} = round(eta_shrinkage_ordered(k0),1);
+ end
+ end
+end
+
+% Separator
+tableCell(end+1,1) = {'<TR>'};
+
+% CORRELATIONS - non-estimated get postfix (FIX).
+for k0=1:length(ALLcorrelationNames),
+ tableCell(end+1,1:2) = {'<TR>' ALLcorrelationNames{k0}};
+ for k=1:length(RESULTS),
+ % Check if parameter in model
+ ix = strmatchIQM(ALLcorrelationNames{k0},RESULTS(k).rawParameterInfo.correlation.names,'exact');
+ if ~isempty(ix),
+ modelVALUES = RESULTS(k).rawParameterInfo.correlation.values(ix);
+ modelRSES = RESULTS(k).rawParameterInfo.correlation.rse(ix);
+ modelESTIMATED = RESULTS(k).rawParameterInfo.correlation.estimated(ix);
+ if modelESTIMATED,
+ value = sprintf('%1.3g',modelVALUES);
+ rse = sprintf('%1.3g%%',modelRSES);
+ valuerse = sprintf('%s %s',postFillCharIQM(value,10,' '),rse);
+ else
+ if ~isnan(modelVALUES),
+ value = sprintf('%1.3g',modelVALUES);
+ rse = sprintf('(FIX)',modelRSES);
+ valuerse = sprintf('%s %s',postFillCharIQM(value,10,' '),rse);
+ else
+ value = NaN;
+ rse = NaN;
+ end
+ end
+ else
+ value = '-';
+ rse = '-';
+ end
+ % Add to table
+ tableCell(end,2+k) = {value};
+ tableCell(end,2+k+1) = {rse};
+ end
+end
+
+% Separator
+tableCell(end+1,1) = {'<TR>'};
+
+% COVARIATES - non-estimated get postfix (FIX).
+for k0=1:length(ALLcovariateNames),
+ tableCell(end+1,1:2) = {'<TR>' ALLcovariateNames{k0}};
+ for k=1:length(RESULTS),
+ % Check if parameter in model
+ ix = strmatchIQM(ALLcovariateNames{k0},RESULTS(k).rawParameterInfo.covariate.names,'exact');
+ if ~isempty(ix),
+ modelVALUES = RESULTS(k).rawParameterInfo.covariate.values(ix);
+ modelRSES = RESULTS(k).rawParameterInfo.covariate.rse(ix);
+ modelESTIMATED = RESULTS(k).rawParameterInfo.covariate.estimated(ix);
+ if modelESTIMATED,
+ value = sprintf('%1.3g',modelVALUES);
+ rse = sprintf('(%1.3g%%)',modelRSES);
+ valuerse = sprintf('%s %s',postFillCharIQM(value,10,' '),postFillCharIQM(rse,14,' '));
+ else
+ if ~isnan(modelVALUES),
+ value = sprintf('%1.3g',modelVALUES);
+ rse = sprintf('(FIX)',modelRSES);
+ valuerse = sprintf('%s %s',postFillCharIQM(value,10,' '),postFillCharIQM(rse,14,' '));
+ else
+ value = NaN;
+ rse = NaN;
+ end
+ end
+ transinfo = getCovCatTransInfo(ALLcovariateNames{k0},RESULTS(k));
+ valuerse = sprintf('%s%s',valuerse,transinfo);
+ else
+ value = '-';
+ rse = '-';
+ transinfo = '-';
+ end
+ % Add to table
+ tableCell(end,2+k) = {value};
+ tableCell(end,2+k+1) = {rse};
+ tableCell(end,2+k+2) = {transinfo};
+ end
+end
+
+% Separator
+tableCell(end+1,1) = {'<TR>'};
+
+% ERROR PARAMETERS - non-estimated get postfix (FIX).
+for k0=1:length(ALLerrorNames),
+ tableCell(end+1,1:2) = {'<TR>' ALLerrorNames{k0}};
+ for k=1:length(RESULTS),
+ % Check if parameter in model
+ ix = strmatchIQM(ALLerrorNames{k0},RESULTS(k).rawParameterInfo.errorParameter.names,'exact');
+ if ~isempty(ix),
+ modelVALUES = RESULTS(k).rawParameterInfo.errorParameter.values(ix);
+ modelRSES = RESULTS(k).rawParameterInfo.errorParameter.rse(ix);
+ modelESTIMATED = RESULTS(k).rawParameterInfo.errorParameter.estimated(ix);
+ if modelESTIMATED,
+ value = sprintf('%1.3g',modelVALUES);
+ rse = sprintf('(%1.3g%%)',modelRSES);
+ valuerse = sprintf('%s %s',postFillCharIQM(value,10,' '),rse);
+ else
+ if ~isnan(modelVALUES),
+ value = sprintf('%1.3g',modelVALUES);
+ rse = sprintf('(FIX)',modelRSES);
+ valuerse = sprintf('%s %s',postFillCharIQM(value,10,' '),rse);
+ else
+ value = NaN;
+ rse = NaN;
+ end
+ end
+ else
+ value = '-';
+ rse = '-';
+ end
+ % Add to table
+ tableCell(end,2+k) = {value};
+ tableCell(end,2+k+1) = {rse};
+ end
+end
+
+% Separator
+tableCell(end+1,1) = {'<TR>'};
+
+% OBJ, AIC, BIC
+tableCell(end+1,1:3) = {'<TR>' 'OBJ' round(RESULTS.OBJ,2)};
+tableCell(end+1,1:3) = {'<TR>' 'AIC' round(RESULTS.AIC,2)};
+tableCell(end+1,1:3) = {'<TR>' 'BIC' round(RESULTS.BIC,2)};
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Display results and save to file
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+% Convert to text and display text
+textDisplay = IQMconvertCellTable2ReportTable(tableCell,'text');
+disp(textDisplay);
+
+% Convert to report text and export to file if filename defined
+IQMconvertCellTable2ReportTable(tableCell,'report',filename);
+
+
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Function to generate transformation and reference value information for
+% covariates ...
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+function [transinfo] = getCovCatTransInfo(ALLcovariateNamesK0,RESULTS)
+COVNAMES = RESULTS.projectHeader.COVNAMES;
+CATNAMES = RESULTS.projectHeader.CATNAMES;
+BETACOVNAMES = RESULTS.projectHeader.BETACOVNAMES;
+BETACOVTRANS = RESULTS.projectHeader.BETACOVTRANS;
+BETACATNAMES = RESULTS.projectHeader.BETACATNAMES;
+BETACATREFERENCE = RESULTS.projectHeader.BETACATREFERENCE;
+
+% Process covariate name to get parameter and covariate name
+name = strrep(ALLcovariateNamesK0,'beta_','');
+name = strrep(name,')','');
+name = strrep(name,'(',',');
+terms = explodePCIQM(name);
+pn = terms{1};
+cn = terms{2};
+% Check if categorical or continuous
+ctypecont = NaN;
+ix = strmatchIQM(cn,COVNAMES,'exact');
+if ~isempty(ix),
+ ctypecont = 1;
+ cc = '';
+else
+ terms = explodePCIQM(cn,'_');
+ cn = terms{1};
+ cc = terms{2};
+ ix = strmatchIQM(cn,CATNAMES,'exact');
+ if ~isempty(ix),
+ ctypecont = 0;
+ end
+end
+if isnan(ctypecont),
+ error('Difficulty to decide covariate type.');
+end
+% Handle continuous covariate
+if ctypecont,
+ % Only need to get the covariate transformation
+ ix = strmatchIQM(ALLcovariateNamesK0,BETACOVNAMES,'exact');
+ trans = BETACOVTRANS{ix};
+ trans = strrep(trans,'cov',cn);
+% transinfo = strrep(sprintf('(%s)',trans),'.','');
+ transinfo = sprintf('%s',trans);
+else
+ % Need to get reference category
+ test = sprintf('beta_%s(%s)',pn,cn);
+ ix = strmatchIQM(test,BETACATNAMES,'exact');
+ trans = BETACATREFERENCE{ix};
+ transinfo = sprintf('Reference group: %s',trans);
+end
+return
diff --git a/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcompareModels.m b/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcompareModels.m
index 9e719ecad2592006ad0867b73a41418492d796aa..46b17f6b9457b57878c12b9556f71b10dad72f1b 100644
--- a/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcompareModels.m
+++ b/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcompareModels.m
@@ -217,16 +217,16 @@ if ~isempty(data),
for k=1:length(catNames),
dataCAT.(catNames{k}) = firstRowData.(catNames{k});
end
- NdataSubjects = length(dataCOV);
+ NdataSubjects = height(dataCOV);
sampleCovariateDosingIX = ceil(NdataSubjects*rand(1,Nsim));
if ~isempty(dataCOV),
- COVvaluesSampled = double(dataCOV(sampleCovariateDosingIX,:));
+ COVvaluesSampled = table2array(dataCOV(sampleCovariateDosingIX,:));
else
COVvaluesSampled = [];
covNames = {};
end
if ~isempty(dataCAT),
- CATvaluesSampled = double(dataCAT(sampleCovariateDosingIX,:));
+ CATvaluesSampled = table2array(dataCAT(sampleCovariateDosingIX,:));
else
CATvaluesSampled = [];
catNames = {};
@@ -376,4 +376,4 @@ else
axis([min(obsTimes) max(obsTimes) get(gca,'YLim')]);
end
-legend(legendText,'Location','NorthWest','Interpreter','none');
\ No newline at end of file
+legend(legendText,'Location','Best','Interpreter','none');
\ No newline at end of file
diff --git a/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcovariateAssessmentUncertainty.m b/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcovariateAssessmentUncertainty.m
index 9c6e1488c3a4345680efc97002694b5219561721..686278a8e3b79626e1cfce4f5a902922af6f1d6d 100644
--- a/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcovariateAssessmentUncertainty.m
+++ b/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcovariateAssessmentUncertainty.m
@@ -542,7 +542,7 @@ for k=1:length(covariateInfo.ParameterName),
tableText{ixRow,ixCol} = sprintf('%1.4g / %1.4g / %1.4g**',quantileIQM(x,0.025),quantileIQM(x,0.5),quantileIQM(x,0.975));
end
end
-tableText(end+1,1:2) = {'<TF>' sprintf('*Parameter value and 95%% uncertainty range for reference subject.\n**Values have been normalized by the median of the non-normalized value (reference subject).')};
+tableText(end+1,1:2) = {'<TF>' sprintf('*Parameter value and 95 percent uncertainty range for reference subject.\n**Values have been normalized by the median of the non-normalized value (reference subject).')};
% Convert to text and display text
textDisplay = IQMconvertCellTable2ReportTable(tableText,'text');
diff --git a/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcreateVPC.m b/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcreateVPC.m
index af69eab857cb20de4d3a57ff95f8b082711e996f..487dc8c67fa97cf8e25e8ca058eed291af5fb667 100644
--- a/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcreateVPC.m
+++ b/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcreateVPC.m
@@ -164,6 +164,21 @@ if ischar(data),
data = IQMloadCSVdataset(data);
end
+% Check if ADDL present and used - in this case: error
+showErrorADDL = 0;
+try
+ ADDL = data.ADDL;
+ % ADDL exists
+ if sum(ADDL>0) ~= 0,
+ showErrorADDL = 1;
+ end
+catch
+ % ADDL does not exist - fine :) - so it does not need to be checked
+end
+if showErrorADDL,
+ error('ADDL column in dataset used to define additional doses. This is not handled in the VPC function.');
+end
+
if ischar(userDosing),
userDosing = IQMdosing(userDosing);
end
diff --git a/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcreateVPCstratified.m b/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcreateVPCstratified.m
index 554f8167b96e56cf554cde43c5a8a92b87168c4e..b5bb3737338d67f65ab4d87bcd09a7e6f33ba776 100644
--- a/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcreateVPCstratified.m
+++ b/IQM Tools/IQMpro/tools/09-NLMEtools/IQMcreateVPCstratified.m
@@ -129,6 +129,21 @@ if ischar(model),
model = IQMmodel(model);
end
+% Check if ADDL present and used - in this case: error
+showErrorADDL = 0;
+try
+ ADDL = data.ADDL;
+ % ADDL exists
+ if sum(ADDL>0) ~= 0,
+ showErrorADDL = 1;
+ end
+catch
+ % ADDL does not exist - fine :) - so it does not need to be checked
+end
+if showErrorADDL,
+ error('ADDL column in dataset used to define additional doses. This is not handled in the VPC function.');
+end
+
% Check GROUP
checkDataColumnsIQM(data,GROUP);
diff --git a/IQM Tools/IQMpro/tools/14-GeneralLinearNLMEmodel/auxiliary/NONMEM/createGeneralLinear_NONMEMprojectIQM.m b/IQM Tools/IQMpro/tools/14-GeneralLinearNLMEmodel/auxiliary/NONMEM/createGeneralLinear_NONMEMprojectIQM.m
index c511f8003c12305f211f8ed89d9894f4d08885b5..ddbc6ee21be2db7c55a764dd32ae7a53abb1e19f 100644
--- a/IQM Tools/IQMpro/tools/14-GeneralLinearNLMEmodel/auxiliary/NONMEM/createGeneralLinear_NONMEMprojectIQM.m
+++ b/IQM Tools/IQMpro/tools/14-GeneralLinearNLMEmodel/auxiliary/NONMEM/createGeneralLinear_NONMEMprojectIQM.m
@@ -246,7 +246,7 @@ try IIVvalues0 = options.IIVvalues0;
try errorModels = options.errorModels; catch, errorModels = ''; end
try errorParam0 = options.errorParam0; catch, errorParam0 = []; end
try covarianceModel = options.covarianceModel; catch, covarianceModel = 'diagonal'; end
-try covariateModel = options.covariateModel; catch, covariateModel = ''; end
+try covariateModel = options.covariateModel; catch, covariateModel = ''; options.covariateModel = ''; end
try covariateModelValues = options.covariateModelValues; catch, covariateModelValues = {}; end
try COVestimate = options.COVestimate; catch, COVestimate = {}; end
@@ -1850,7 +1850,7 @@ COVARIATENAMES_info = sprintf('; COVARIATENAMES = ''%s''\r\n',x(1:end-1));
PROJECT_INFO_TEXT = sprintf('%s%s',PROJECT_INFO_TEXT,COVARIATENAMES_info);
% COVARIATESUSED
-COVARIATESUSED = setdiff(explodePCIQM(strrep(strrep(covariateModel,'{',''),'}','')),parameterNames);
+COVARIATESUSED = setdiff(explodePCIQM(strrep(strrep(options.covariateModel,'{',''),'}','')),parameterNames);
x = sprintf('%s,',COVARIATESUSED{:});
COVARIATESUSED_info = sprintf('; COVARIATESUSED = ''%s''\r\n',x(1:end-1));
PROJECT_INFO_TEXT = sprintf('%s%s',PROJECT_INFO_TEXT,COVARIATESUSED_info);