% figure 4 (revision)
%% A

d = obs_all(goodGenes, :);
n = length(find(d(:, 2) > d(:, 1)))*100./length(goodGenes);
disp(['between hybrid alleles > between species: ', num2str(n), '%']);

comparativeTh = 0.2;
controlTh = 0.5;

thSpecies = comparativeTh;
thHybrid = controlTh;
thExpHyb = controlTh;
thCis = comparativeTh;
% relativeToDiagTh = -.8;

n1 = length(find(obs_all(goodGenes, 1) < thSpecies));
transGenes = find(obs_all(goodGenes, 1) < thSpecies & obs_all(goodGenes, 2) > thHybrid);
disp(['Rc species < ', num2str(thSpecies), ' = ', num2str(n1)]);
disp(['trans = ', num2str(length(transGenes))]);

cisGenes = find(obs_all(goodGenes, 2) < thCis & expH(goodGenes) > thExpHyb);
% relativeToDiag = obs_all(goodGenes, 2) - expH(goodGenes);
% cisGenes = find(obs_all(goodGenes, 2) < thCis & expH(goodGenes) > thExpHyb & relativeToDiag < relativeToDiagTh);

exampleGeneNames = {'VTC3', 'LSO2', 'CHA1'};
% exampleGenes = {'VTC3', 'CHA1', 'LSO2', 'PGM1'};
exampleGenes = ismember_smart(allGenes(intint(goodGenes)), exampleGeneNames);

fs = 14;
figure;
subplot(1,2,1);
x = obs_all(goodGenes, 1); 
y = obs_all(goodGenes, 2);
dscatter(x, y);
hold on; 
cm = cbrewer2('YlGnBu'); cm = cm(end:-1:1, :);
colormap(cm);
axis square; axis([-1 1 -1 1]); plot(xlim, ylim, 'k'); 
set(gcf, 'color', 'w', 'position', [263   302   972   577]); box on; grid on;
xlabel('Regulatory similarity, species'); ylabel('Regulatory similarity, hybrid');
x1 = [-.8 thSpecies thSpecies -.8];
y1 = [thHybrid thHybrid .95 .95];
patch(x1, y1, 'w', 'EdgeColor', 'r', 'FaceColor', 'none');
set(gca, 'fontsize', fs);
scatter(x(exampleGenes(1)), y(exampleGenes(1)), [], rgb('LimeGreen'), 'filled');
scatter(x(exampleGenes(2)), y(exampleGenes(2)), [], rgb('Orange'), 'filled');
scatter(x(exampleGenes(3)), y(exampleGenes(3)), [], rgb('Pink'), 'filled');
text(x(exampleGenes), y(exampleGenes), exampleGeneNames);
expH = nanmedian(exp_all(:, 3:4), 2);

x = expH;
y = obs_all(:, 2);
x = x(goodGenes); y = y(goodGenes);
subplot(1,2,2); 
dscatter(x, y); hold on;
colormap(cm)
plot([0 1], [0 1], 'k');
xlabel('dataset control');
ylabel('Regulatory similarity, hybrid');
axis square; box on; grid on; set(gcf, 'color', 'w');
x1 = [thExpHyb 1 1 thExpHyb];
y1 = [-1 -1 thCis thCis];
patch(x1, y1, 'w', 'EdgeColor', 'r', 'FaceColor', 'none');
scatter(x(exampleGenes(1)), y(exampleGenes(1)), [], rgb('LimeGreen'), 'filled');
scatter(x(exampleGenes(2)), y(exampleGenes(2)), [], rgb('Orange'), 'filled');
scatter(x(exampleGenes(3)), y(exampleGenes(3)), [], rgb('Pink'), 'filled');
text(x(exampleGenes), y(exampleGenes), exampleGeneNames);
set(gca, 'fontsize', fs);

g = intint(goodGenes(cisGenes));
a = allGenes(g);
b = GP.gene_infoR64.desc(g);
t = table(a,b, 'VariableName', {'name', 'description'});
disp(['# trans genes: ', num2str(length(transGenes)), ', # cis genes : ', num2str(length(cisGenes))]);

Rc_cis = intint(goodGenes(cisGenes));
Rc_trans = intint(goodGenes(transGenes));
Rc_genes = struct;
Rc_genes.Rc_cis = Rc_cis;
Rc_genes.Rc_trans = Rc_trans;
Rc_genes.comparativeTh = comparativeTh;
Rc_genes.contorlTh = controlTh;
save('Rc_genes.mat', '-struct', 'Rc_genes');

%% D / F

fs = 12;
cisTransInAbs = medianEffect(:, 2:3);
cisTransInAbs = cisTransInAbs(intint, :);
cisTransInAbs = cisTransInAbs(goodGenes, :);
cisTransInAbs = abs(cisTransInAbs);

g = transGenes;
currLabel = 'trans';
figure3_panelD;
currLabel = 'cis';
g = cisGenes;
figure3_panelD;