nav[aria-label="Primary Navigation"] { padding: 0; & ul { list-style: none; width: 100%; display: flex; flex-direction: row; justify-content: start; align-items: start; gap: 30px; padding: 0; & li { margin: 0; } & ul li { list-style: none; } } }

Extract factor range from bivariate graph

pauldeen

How do I efficiently extract the crossing points of this curve and its 95% CI with an arbitrary Y value.

dt = Open( "$sample_data/Growth.jmp" );
biv = Bivariate(
	Y( :age ),
	X( :ratio ),
	Fit Line( {Line Color( {204, 121, 41} )}, save studentized residuals ),
	SendToReport(
		Dispatch(
			{},
			"Bivar Plot",
			FrameBox,
			{Grid Line Order( 3 ), Reference Line Order( 4 )}
		)
	)
);
biv << (Curve[1] << Save Studentized Residuals);

biv2 = Bivariate(
	Y( :Studentized Residuals age ),
	X( :ratio ),
	Fit Polynomial( 2, {Confid Shaded Fit( 1 ), Line Color( {212, 73, 88} )} ),
	SendToReport(
		Dispatch(
			{},
			"ratio",
			ScaleBox,
			{Min( 0.2796875 ), Max( 1.2666259765625 ), Inc( 0.2 ), Minor Ticks( 1 )}
		),
		Dispatch(
			{},
			"Studentized Residuals age",
			ScaleBox,
			{Add Ref Line( 1.96, "Solid", "Black", "", 1 ),
			Add Ref Line( -1.96, "Solid", "Black", "", 1 )}
		)
	)
);

So I would like to get the Ratio values where the curve and its 95%CI crosses +/- 1.96 (so 6 values in total). Are there any better solutions than the brute force I can think of: dumping the predicteds and the mean confidence limit formula into a table and adding rows with values to the table to see when it crosses the limits?

Accepted answers

All comments

I don't know of a direct method but here are a few ideas. Some would be a bit of work to script and probably I suspect none would be much faster for small problems.

Save a formula for each of those lines to the data table by pulling the parameters and standard errors from the table in the bivariate report and then creating a column with the appropriate formula. Then use the optimizer built into the profiler to solve for each.
Save relatively few points using the brute force method and then fit a new curve to them using fit model, and then save the prediction formula and solve with the optimizer.
Fit using fit model instead where you should be able to extract the prediction formula directly, but maybe not the CI.
Use the brute force method and interpolate between the closest few points, either linear between the last two, or for this problem quadratic between the closest three.

Mark_Bailey

You might be able to use the Minimize() function and the expression for the confidence boundary. See the Help for optimization methods.

pauldeen

@Mark_Bailey thanks, I tried giving that a shot but the syntax is unclear to me.

This is as far as I got, but now I need to figure out how to feed it it the right optimization constraints:

Names Default To Here( 1 );

dt = Open( "$sample_data/Growth.jmp" );
biv = Bivariate(
	Y( :age ),
	X( :ratio ),
	Fit Line( {Line Color( {204, 121, 41} )}, save studentized residuals ),
	SendToReport(
		Dispatch(
			{},
			"Bivar Plot",
			FrameBox,
			{Grid Line Order( 3 ), Reference Line Order( 4 )}
		)
	)
);
biv << (Curve[1] << Save Studentized Residuals);
biv << close window();

biv2 = Bivariate(
	Y( :Studentized Residuals age ),
	X( :ratio ),
	Fit Polynomial( 2, {Confid Shaded Fit( 1 ), Line Color( {212, 73, 88} )} ),
	SendToReport(
		Dispatch(
			{},
			"ratio",
			ScaleBox,
			{Min( 0.2796875 ), Max( 1.2666259765625 ), Inc( 0.2 ), Minor Ticks( 1 )}
		),
		Dispatch(
			{},
			"Studentized Residuals age",
			ScaleBox,
			{Add Ref Line( 1.96, "Solid", "Black", "", 1 ),
			Add Ref Line( -1.96, "Solid", "Black", "", 1 )}
		)
	)
);
biv << (Curve[1] << Mean Confidence Limit Formula);
f =column(dt, n cols(dt)) << Get Formula;
fs = substitute(name expr(f), expr(:ratio), expr(x1));
/*A = [1 5];
b = [-1.96, 1.96];
minFun = Constrained Maximize(
	f,
	{x1( 0, 5 )},
	<<lessthanEQ( {A, b} )/*and/or <<GreaterThanEQ({A,b}) and/or <<EqualTo({A,b}),
	<<StartingValues( [1, .5] )
);
Eval List( {x1, minFun} );*/<BR />
minimize(name expr(f), {})

The last part doesn't work. Thanks for your help!

Mark_Bailey

I wonder about another approach. What if you use the Prediction Profiler and a Match Target desirability function with your target response? Then you can just pick off the predicted mean and the 95% confidence interval estimate of the mean.

I might have strayed from your original goal...

Mark_Bailey

Your code was a bit messy. I cleaned it up and finished the example, though I might still be going after the wrong goal.

Names Default To Here( 1 );

// example data set
dt = Open( "$sample_data/Growth.jmp" );

// fit first-order polynomial model to response
biv = Bivariate(
	Y( :age ),
	X( :ratio ),
	Fit Line
);

// save residuals for next step
biv << (Curve[1] << Save Studentized Residuals);
biv << Close Window();

// fit second-order polynomial model to residuals
biv = Bivariate(
	Y( :Studentized Residuals age ),
	X( :ratio ),
	Fit Polynomial( 2, {Confid Shaded Fit( 1 ), Line Color( {212, 73, 88} )} )
);

// save formulas for lower / upper 95% confidence bound
biv << (Curve[1] << Mean Confidence Limit Formula);

// get upper formula
f = Column( dt, N Cols( dt ) ) << Get Formula;

// replace predictor name with x
fs = Substitute( Name Expr( f ), Expr( :ratio ), Expr( x ) );

// initialize x, could get x that predicts mean residual = 0
x = 0.5;

// minimize expression that is Abs( residual ), or mean = 0
fa = Insert( Expr( Abs() ), Name Expr( fs ) );

Minimize( fa, { x } );
Show( x );

pauldeen

Yeah the problem with that one is that there are usually multiple crossings with the target value (one on the left and one on the right) and with the desirability maximizer you randomly get one of the two.

I tried Inverse predictions under fit model but that does not support quadratic terms.

pauldeen

Yeah, exactly this works for finding the local minimum but it still is hard to get to the next crossing with 1.96 from there. You would still need to walk up the curve untill values exceed 1.96.

To temperoraly solve this, I built function for it that just takes little steps and goes out in both directions untill the value exceeds the limits (or the max number of steps) but I had expected there to be a solver function in JMP that would implement a better algorithm (like what is built into the non-linear platform). I guess for now I have what I need but I'm always open for nicer ways to do it.