% Given the following variables: 
% 'vqm' is a column vector that holds one objective model scores for the superset. 
% 'vqm_sign' = 1 or -1 and gives the direction of 'vqm' with respect to the subjective scale. 
% Higher values of vqm for higher subjective scores means that 'vqm_sign' = 1. 
% 'mos' is a column vector that holds the corresponding MOS values for the superset. 
% 'order' specifies the order of the monotonic polynomial fit (e.g., order=3). 
% 
% Following code implements monotonic polynomial fitting using the MATLAB Optimization 
% Toolbox routine lsqlin. 
% 
 
num_comb = length(vqm); % total number of clips in the superset 
 
% Create x and dx arrays. The dx slope array (holds the derivatives of mos with respect 
% to vqm), the vqm_sign specifies the direction of the slope that must not change over the 
% vqm range. 
x = ones(num_comb,1); 
dx = zeros(num_comb,1); 
for col = 1:order
    x = [x vqm.^col];
    dx = [dx col*vqm.^(col-1)];
end
 
% The lsqlin routine uses <= inequalities. Thus, if vqm_sign is -1 (negative slope), 
% we are correct but if vqm_sign is +1 (positive slope), we must multiply each side by -1. 
if (vqm_sign == 1)
    dx = -1*dx;
end
 
% Perform the fit and organize to what would have been output by polyfit 
fit = lsqlin(x,mos,dx,zeros(num_comb,1)); 
fit = flipud(fit)'; 
 
% Use polyval to find the fitted vqm values (e.g., vqm_hat) 
vqm_hat = polyval(fit,vqm); 