@dharr yes, sure about the threshold line. Recap:

1. I want to see where A>B, with A=sigma_d*(2*n-1)/(2*sqrt(4*n-3)) and B=sigma_dr/sqrt(n). Given the positivity assumptions, this reduces to showing where X=4*sigma_d^2*n^3-4*sigma_d^2*n^2+(sigma_d^2-16*sigma_dr^2)*n+12*sigma_dr^2>0.
    
2. Divide by sigma_d^2 and let y=sigma_dr/sigma_d and X reduces to XX=4*n^3-4*n^2+(1-16*y^2)*n+12*y^2.
    
3. Since n>=2, y>=3/sqrt(10)~0.95. Solve XX=0 for y: y*=((2*n-1)/2)*sqrt(n/4*n-3).
    
4. Conclusion:
   - Unallowed region: y<0.95
   - Allowed regions: for y>=0.95, if y<y* then XX>0 otherwise XX<0 (from your 3D plot, am I right?)

Question: what do I conclude about A and B? Not that A>B if 0.95 <= sigma_dr/sigma_d < y* right? Don't I need to scale to go back to my original X (and hence to A and B) from XX? This is what I meant by "going back", sorry for the misunderstanding. So we know where XX>0 now, but I want to know in which region X>0.

@dharr is this correct then and do I need interpolation at all to find the threshold line? If yes, how do I go back to my parameters plot for x from the plot for x^2

Finally, is this consistent with the CellDecomposition for n>0

Worksheet: regions_quartic_equation_in_n_MaPal_more.mw

From CellDecomposition, this cubic in n seems to have only one root for n>=2...Now If I choose sigma__dr=sigma__d it's easy to show that for n>=2 the curve is positive for n>2.1 only (that is, for n=3,4,5...). But can I show it more generally as you did for the original quartic? quartic_equation_in_n_MaPal.mw 

@Kitonum the two real roots would be these quartic_parametric.mw? 

What if n is constrained to be >1 and integer? no roots seem to exist. Is this trivial? Is my equation always positive and, if so, how do I show this analytically?

@dharr sorry I forgot to specify that.

n>=2 and must be an integer. So n=2,3,4...I edited the question accordingly. 

Here I apply CellDecomposition, but it doesn't work if I choose n>=2 instead of n>0 as constraint:
quartic_parametric.mw

Is it because for n>=2 my equation is always positive (no roots)? If so, how could I prove this analytically?  

@dharr @C_R, I was wondering if you had any update on this question

@dharr good that you were thinking something similar.

Assuming a rearrangement into a form containing only L (and X) exists, does it mean that theoretically from such L we can go back and find l[i] from each eqi for i from 1 to any n (and not just n=3)?

I thought this because of the similarity of eqi as well as the L they all share...almost like if any eqi can be rewritten as a x[i]*K where K is the same across eqi...  

@C_R I would avoid numeric...

Question (also for @dharr): Can I rewrite l[1], l[2], l[3], and l__sq.w__sq as some functions of L, and then solve {eqq1, eqq2, eqq3} for L?
See attached: MaPal_attempt_dot_products_for_C_R_reply.mw. Note that L is defined as l.w__sq and, as such, it appears in all three equations!

The idea is to keep everything as implicit as possible and not working with components, so that I can perhaps obtain a closed form solution for lambda (or lambda^T once I evaluate it as a vector at the very end) rather than for lambda[1], lambda[2], and lambda[3]. Do you see my point?

@dharr thanks.

I tried to set up the problem more in matrix form, but it's not helpful so far: attempt_dot_products_for_dharr.mw.

The numerator of my equations is a covariance, while the denominator is a variance. I computed these by hand. I expressed as matrices or vectors all that I could express as such. In particular, R is a matrix such that R[i,i]=r[i] and R[i,j]=0 (i<>j), while Delta is a matrix such that Delta[i,i]=delta[i] and Delta[i,j]=0 (i<>j). r[i] and delta[i] are defined in the worksheet.

About r vector. It is made of ratios of standard deviations and is derived from somewhere else. Why would its origin matter? Cannot it be kept implicit (e.g., sq=sqrt(A/B)) through the end since its a scalar?

Overall, I understand that you recommend working with a consistent matrix formulation from start to finish rather than a hybrid. But I am not sure if a pure matrix formulation even exists. What exactly do you think is missing from attempt_dot_products_for_dharr.mw.? If you are specific, I can see if I can come up with a more matrixy formulation for the elements you identify ...

About avoiding working with components. I haven't tried this yet, but I see your point.

@C_R after ~1.5h still evaluating...

Worksheet: attempt_dot_products_for_C_R.mw

@C_R thanks for your reply.

I tried to address your points, but I failed. Please check: Attempt_dot_products.mw

@dharr thank you for clarifying what you meant about the calculation errors.

Also thank you for re-doing the calculations! Yes, now it's definitely correct.

@dharr thanks for the details on that limit. About the signs: the limit for lambda must be positive, but no problem if the resulting limit for X ends up being negative for some Delta, d, and w. That's fine.

So back to the original question: given that the X are indeed the same for the two approaches, why the evaluations for fixed n look different in the two scripts?

Here's what I find: splitting_into_two.mw
W_XiXj_only captures the exclusive contribution of X[i] and X[j] on W. I compute this by setting X__r=0. W_Xr_only captures the exclusive contribution of X__r on W. I compute this by setting X[i]=X[j]=0. I was expecting  W_XiXj_only+W_Xr_only=W but this is not what I find. Mathematically, is it because X[i]*X__r (and X[j]*X__r) terms are 0 in both exclusive cases so when I sum them up I am neglecting these "interaction terms"? If so, why is this evident when I leave limit lambda undefined until the very end (your approach) but not when I define it at the onset?

n=2 seems to be an exception. In this case there are no "interaction terms" anyway.

@dharr here: same.mw

(same X).

Sorry for not including this earlier.

@dharr thank you for your help. Why evaluating your simple result for n={3,4,5,6} do not coincide with the previous evaluations? (for n=2 they do)

Here the evaluation of your result after defining the (gamma to infinity limit of) lambda: generic_n-form_of_infinity_limit2_MaPal.mw
(the gamma to infinity limits of X[i] and X[j] are automatically taken care of)

Here the previous evaluation where the (gamma to infinity limits of) lambda, X[i], and X[j] are already embedded in the definition of W before computing its limit: n_from_2_to_6_-_not_coincide.mw (but with your d-d_i and d-d_j trick)

Mathematically (and sorry if this is not a Maple question), shouldn't the two evaluations be equivalent? Or are they equivalent only if limit(W,gamma=infinity) is continuous at limit(X[i/j],gamma=infinity)?

Overall, I am more confident in your approach where I leave X[i] and X[j] as functions of gamma first, then compute the limit of the whole W, and only at the very end define the lambda. However, some things I can't explain. For example, I was not expecting delta_r to multiply other deltas (d) in the last term of your result...