Se solicita ayuda con scripts: creación de un simulador para errores de tipo 1 en pruebas de equivalenciaHola a todos,
Estoy intentando crear un script que genere un conjunto de datos aleatorio (n = 3 mediciones de referencia y de prueba) y luego ejecute una simulación con un número determinado (1000 o el que sea necesario). Para cada conjunto, ejecutaré una prueba de equivalencia y me gustaría que el resultado sea un informe/tabla con la frecuencia de los eventos donde se demostró la equivalencia. El objetivo de esta herramienta será proporcionar una probabilidad de errores de tipo 1 dado un cierto margen de equivalencia.
He encontrado el siguiente código, pero no funciona y estoy prácticamente atascado. Agradecería cualquier ayuda que el grupo pueda brindarme.
¡Gracias!
Chris
// Parameters
numSimulations = 10;
nPerGroup = 3;
refMean = 99.1;
testMean = 98.8;
sd = 0.1;
equivMargin = 0.3;
type1Errors = 0;
// Critical t-value for 90% CI with df = 2n - 2
df = 2 * nPerGroup - 2;
tCrit = StudentTQuantile(0.950, df);
// Create summary table
summaryTable = New Table("Simulation Summary",
New Column("Simulation", Numeric),
New Column("CI_Lower", Numeric),
New Column("CI_Upper", Numeric),
New Column("Equivalence_Concluded", Character),
New Column("Mean_Diff", Numeric)
);
// Loop through simulations
For(sim = 1, sim <= numSimulations, sim++,
// Create new data table
dt = New Table("Sim" || Char(sim),
Add Rows(nPerGroup * 2),
New Column("Group", Character),
New Column("Value", Numeric)
);
// Simulate Reference group
For(i = 1, i <= nPerGroup, i++,
dt:Group[i] = "Reference";
dt:Value[i] = Random Normal(refMean, sd);
);
// Simulate Test group
For(i = nPerGroup + 1, i <= nPerGroup * 2, i++,
dt:Group[i] = "Test";
dt:Value[i] = Random Normal(testMean, sd);
);
// Extract values
refVals = dt << Get As Matrix( :Value )[1::nPerGroup];
testVals = dt << Get As Matrix( :Value )[nPerGroup+1::nPerGroup*2];
meanDiff = Mean(testVals) - Mean(refVals);
// Sample variances
varRef = Variance(refVals);
varTest = Variance(testVals);
// Pooled variance
pooledVar = ((nPerGroup - 1) * varRef + (nPerGroup - 1) * varTest) / df;
pooledSD = Sqrt(pooledVar);
// Standard error of difference
seDiff = pooledSD * Sqrt(2 / nPerGroup);
// 90% CI
ciLower = meanDiff - tCrit * seDiff;
ciUpper = meanDiff + tCrit * seDiff;
// Check equivalence
eqResult = "No";
If(ciLower > -equivMargin & ciUpper < equivMargin,
type1Errors++;
eqResult = "Yes";
);
// Log results
summaryTable << Add Rows(1);
summaryTable:Simulation[sim] = sim;
summaryTable:CI_Lower[sim] = ciLower;
summaryTable:CI_Upper[sim] = ciUpper;
summaryTable:Equivalence_Concluded[sim] = eqResult;
summaryTable:Mean_Diff[sim] = meanDiff;
// Close the table to avoid clutter
dt << Close Window;
);
// Show results
New Window("Type I Error Summary",
Text Box("Total Simulations: " || Char(numSimulations)),
Text Box("Type I Errors: " || Char(type1Errors)),
Text Box("Estimated Type I Error Rate: " || Char(Round(type1Errors / numSimulations * 100, 2)) || "%")
);
Automatización y scriptingRe: Se solicita ayuda con scripts: creación de un simulador para errores de tipo 1 en pruebas de equivalencia¿En qué parte estás atascado?
Re: Se solicita ayuda con scripts: creación de un simulador para errores de tipo 1 en pruebas de equivalenciaH. Jarmo,
Actualmente, el código falla al intentar extraer los valores T críticos del informe de equivalencia. Sin embargo, existen otros problemas. Lo he tenido en un estado en el que ejecuta la simulación y las pruebas de equivalencia, pero no puede extraer los valores del informe.
También intenté aplicar algo de esto a la IA, pero creo que algunas de las convenciones de sintaxis y codificación pueden ser incorrectas.
Chris
Re: Se solicita ayuda con scripts: creación de un simulador para errores de tipo 1 en pruebas de equivalenciaSugeriría estas ediciones:
No existe la función StudentTQuantile(). Debería ser t Quantile().
Variance() no es una función integrada. Defínela en su script (o simplemente eleve al cuadrado Std Dev()).
refVals = dt << Get As Matrix( :Value )[1::nPerGroup] debería ser más como refVals = dt:Value[Where( :Group == "Reference" )];
Lo mismo para testVals =
Establezca la visibilidad de dt en "Privado" para acelerar considerablemente las cosas.
Variance = Function( {Values},
Return( Std Dev( Values ) ^ 2 )
);
// Parameters
numSimulations = 1000;
nPerGroup = 3;
refMean = 99.1;
testMean = 98.8;
sd = 0.1;
equivMargin = 0.3;
type1Errors = 0;
// Critical t-value for 90% CI with df = 2n - 2
df = 2 * nPerGroup - 2;
tCrit = t Quantile( 0.950, df );
// Create summary table
summaryTable = New Table( "Simulation Summary",
New Column( "Simulation", Numeric ),
New Column( "CI_Lower", Numeric ),
New Column( "CI_Upper", Numeric ),
New Column( "Equivalence_Concluded", Character ),
New Column( "Mean_Diff", Numeric )
);
// Loop through simulations
For( sim = 1, sim <= numSimulations, sim++,
// Create new data table
dt = New Table( "Sim" || Char( sim ), "Private", Add Rows( nPerGroup * 2 ), New Column( "Group", Character ), New Column( "Value", Numeric ) );
// Simulate Reference group
For( i = 1, i <= nPerGroup, i++,
dt:Group[i] = "Reference";
dt:Value[i] = Random Normal( refMean, sd );
);
// Simulate Test group
For( i = nPerGroup + 1, i <= nPerGroup * 2, i++,
dt:Group[i] = "Test";
dt:Value[i] = Random Normal( testMean, sd );
);
// Extract values
refVals = dt:Value[Where( :Group == "Reference" )];
testVals = dt:Value[Where( :Group == "Test" )];
meanDiff = Mean( testVals ) - Mean( refVals );
// Sample variances
varRef = Variance( refVals );
varTest = Variance( testVals );
// Pooled variance
pooledVar = ((nPerGroup - 1) * varRef + (nPerGroup - 1) * varTest) / df;
pooledSD = Sqrt( pooledVar );
// Standard error of difference
seDiff = pooledSD * Sqrt( 2 / nPerGroup );
// 90% CI
ciLower = meanDiff - tCrit * seDiff;
ciUpper = meanDiff + tCrit * seDiff;
// Check equivalence
eqResult = "No";
If( ciLower > -equivMargin & ciUpper < equivMargin,
type1Errors++;
eqResult = "Yes";
);
// Log results
summaryTable << Add Rows( 1 );
summaryTable:Simulation[sim] = sim;
summaryTable:CI_Lower[sim] = ciLower;
summaryTable:CI_Upper[sim] = ciUpper;
summaryTable:Equivalence_Concluded[sim] = eqResult;
summaryTable:Mean_Diff[sim] = meanDiff;
// Close the table to avoid clutter
dt << Close Window;
);
// Show results
New Window( "Type I Error Summary",
Text Box( "Total Simulations: " || Char( numSimulations ) ),
Text Box( "Type I Errors: " || Char( type1Errors ) ),
Text Box( "Estimated Type I Error Rate: " || Char( Round( type1Errors / numSimulations * 100, 2 ) ) || "%" )
);
Re: Se solicita ayuda con scripts: creación de un simulador para errores de tipo 1 en pruebas de equivalenciaJarmo,
Gracias por la sugerencia y la explicación. ¡Funciona de maravilla!
Voy a trabajar con el código para intentar internalizar cómo funciona todo... ¡Agradezco tu ayuda!
Chris
Re: Se solicita ayuda con scripts: creación de un simulador para errores de tipo 1 en pruebas de equivalenciaEl código original parece bastante probable que haya sido creado por LLM (debido a los errores obvios que cometen si conoces algo de JSL).
Ya se ha proporcionado una versión funcional, pero así es como yo lo abordaría. Esta es una versión un poco más compleja que utiliza matrices (generalmente son más rápidas) y no he verificado los cálculos.
Names Default To Here(1);
numSimulations = 1000;
nPerGroup = 3;
refMean = 99.1;
testMean = 98.8;
sd = 0.1;
equivMargin = 0.3;
type1Errors = 0;
df = 2 * nPerGroup - 2;
tCrit = t Quantile(0.950, df);
dt = New Table("Simulation Summary",
New Column("Simulation", Numeric),
New Column("CI_Lower", Numeric),
New Column("CI_Upper", Numeric),
New Column("Equivalence_Concluded", Character),
New Column("Mean_Diff", Numeric),
Invisible
);
For(sim = 1, sim <= numSimulations, sim++,
m_ref = J(nperGroup, 1, Random Normal(refMean, sd));
m_test = J(nperGroup, 1, Random Normal(testMean, sd));
avgdiff = Mean(m_test) - Mean(m_ref);
var_ref = Std Dev(m_ref)^2;
var_test = Std Dev(m_test)^2;
pooledVar = ((nPerGroup - 1) * var_ref + (nPerGroup - 1) * var_test) / df;
pooledSD = Sqrt(pooledVar);
seDiff = pooledSD * Sqrt(2 / nPerGroup);
ciLower = avgdiff - tCrit * seDiff;
ciUpper = avgdiff + tCrit * seDiff;
eqResult = "No";
If(ciLower > -equivMargin & ciUpper < equivMargin,
type1Errors++;
eqResult = "Yes";
);
dt << Add Rows(
{Simulation = sim, CI_Lower = ciLower, CI_Upper = ciUpper, Equivalence_Concluded = eqResult, Mean_Diff = avgdiff}
);
);
dt << Show Window(1);
nw = New Window("Type I Error Summary",
Text Box("Total Simulations: " || Char(numSimulations)),
Text Box("Type I Errors: " || Char(type1Errors)),
Text Box("Estimated Type I Error Rate: " || Char(Round(type1Errors / numSimulations * 100, 2)) || "%")
);