image thumbnail
from Linear regression for multiple category data with different slopes for each by Soumya Banerjee
The function does OLS and RMA regression on data that is categorized into two different categories

regress_plot_twocategories_allslope

Contents

function regress_plot_twocategories_allslope(input_x,input_y,input_x1,input_y1,input_x2,input_y2,desc,x_label,y_label)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Name - regress_plot_twocategories_allslope
% Creation Date - 18th Oct 2011
% Author: Soumya Banerjee
% Website: www.cs.unm.edu/~soumya
%
% Acknowledgements: Dr. Melanie Moses
%
% Description: Takes six column vectors, a description, x label and y label and
% plots the data and outputs all the statistics (r2, OLS slope, RMA slope, 95% CI
% intervals)
% Same syntax as plot(x,y) : x-axis data 1st followed by y-axis data
% Plots the slopes of both the categories as well as combined
%
% input_x:  combined column vector x
% input_y:  combined column vector y
% input_x1: category 1 column vector x
% input_y1: category 1 column vector y
% input_x2: category 2 column vector x
% input_y2: category 2 column vector y
% desc:     description of plot (goes in the title of the plot)
% x_label:  label for x-axis of plot
% y_label:  label for y-axis of plot
%
% Example usage:
% regress_plot_twocategories_allslope([1 1 1 1 2 2 2 2 2]',[1 1.1 1.2 1.3 2.2 2.1 2.3 2 2.4]',...
% [1 1 1 1]', [1 1.1 1.2 1.3]', [2 2 2 2 2]',[2.2 2.1 2.3 2 2.4]',...
% 'analysis of all individuals mcmc burn-in thinning',...
% 'log_1_0(M)','log_1_0(delta)')
%
% License - BSD
%
% Change History -
%                   18th Oct 2011 - Creation by Soumya Banerjee
%                   5th  Jan 2013 - Modified by Soumya Banerjee
%                                       Plots slopes for both categories
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

combined

msapop_ones = [ones(size(input_x)), input_x];
[b1 bint r rint stats] = regress(input_y,msapop_ones);
display('RMA Slope')
rmaslope = b1(2,1) / stats(1,1) ^ 0.5
display('bint')
bint
display('95% CI on intercept')
bint(1,:)
display('95% CI on slope')
bint(2,:)
display('r2')
stats(1,1)
display('OLS Slope')
b1(2,1)
olsslope = b1(2,1);
display('p value')
stats(1,3)
display('b1')
b1
Error using regress_plot_twocategories_allslope (line 43)
Not enough input arguments.

first category

msapop_ones_1 = [ones(size(input_x1)), input_x1];
[b1_1 bint_1 r_1 rint_1 stats_1] = regress(input_y1,msapop_ones_1);
display('RMA Slope category 1')
rmaslope_1 = b1_1(2,1) / stats_1(1,1) ^ 0.5
display('bint category 1')
bint_1
display('95% CI on intercept category 1')
bint_1(1,:)
display('95% CI on slope category 1')
bint_1(2,:)
display('r2 category 1')
stats_1(1,1)
display('OLS Slope category 1')
b1_1(2,1)
olsslope_1 = b1_1(2,1);
display('p value category 1')
stats_1(1,3)
display('b1 category 1')
b1_1

second category

msapop_ones_2 = [ones(size(input_x2)), input_x2];
[b1_2 bint_2 r_2 rint_2 stats_2] = regress(input_y2,msapop_ones_2);
display('RMA Slope category 2')
rmaslope_2 = b1_2(2,1) / stats_2(1,1) ^ 0.5
display('bint category 2')
bint_2
display('95% CI on intercept category 2')
bint_2(1,:)
display('95% CI on slope category 2')
bint_2(2,:)
display('r2 category 2')
stats_2(1,1)
display('OLS Slope category 2')
b1_2(2,1)
olsslope_2 = b1_2(2,1);
display('p value category 2')
stats_2(1,3)
display('b1 category 2')
b1_2

plot

figure
plot(input_x1,input_y1,'ks','MarkerEdgeColor','k','MarkerFaceColor','k','MarkerSize',10)
xlabel(x_label,'fontsize',23);
ylabel(y_label,'fontsize',23);
title(strcat(strcat(strcat(strcat(desc,strcat(strcat(strcat(strcat('r^2 value = ',num2str(stats(1,1))), strcat('OLSslope = ',num2str(olsslope))), 'RMAslope = '),num2str(rmaslope)))),'p value = '),num2str(stats(1,3))),'fontsize',18);
hold on
plot(input_x2,input_y2,'ro','MarkerEdgeColor','k','MarkerFaceColor','r','MarkerSize',10)
%  axis([-2 0.55 0.25 0.75])
% axis([1 11 1 17])
plot(input_x,olsslope .* input_x + b1(1,1),'-b');
hold on
% plot([-2 -1 0 1 2 3],[5 5 5 5 5 5],'-r')
% plot(plot_x,rmaslope * plot_x + 1,'-g');
plot(input_x1,olsslope_1 .* input_x1 + b1_1(1,1),'-k');
hold on
plot(input_x2,olsslope_2 .* input_x2 + b1_2(1,1),'-r');
hold on
legend('Category 1','Category 2','Combined','Location','NorthEast')

hold off

display('b1')
b1
end

Contact us