xlabel name legends font change achsenbeschriftung matlab r statistics

name - title image matlab



Prueba de independencia de Matlab (2)

Para 1,000,000 de observaciones, observé un evento discreto, X, 3 veces para el grupo de control y 10 veces para el grupo de prueba.

Necesito realizar una prueba de independencia de Chi cuadrado en Matlab. Así es como lo harías en r:

m <- rbind(c(3, 1000000-3), c(10, 1000000-10)) # [,1] [,2] # [1,] 3 999997 # [2,] 10 999990 chisq.test(m)

La función r devuelve chi-squared = 2.7692, df = 1, p-value = 0.0961.

¿Qué función de Matlab debo usar o crear para hacer esto?


Aquí está mi propia implementación que uso:

function [hNull pValue X2] = ChiSquareTest(o, alpha) %# CHISQUARETEST Pearson''s Chi-Square test of independence %# %# @param o The Contignecy Table of the joint frequencies %# of the two events (attributes) %# @param alpha Significance level for the test %# @return hNull hNull = 1: null hypothesis accepted (independent) %# hNull = 0: null hypothesis rejected (dependent) %# @return pValue The p-value of the test (the prob of obtaining %# the observed frequencies under hNull) %# @return X2 The value for the chi square statistic %# %# o: observed frequency %# e: expected frequency %# dof: degree of freedom [r c] = size(o); dof = (r-1)*(c-1); %# e = (count(A=ai)*count(B=bi)) / N e = sum(o,2)*sum(o,1) / sum(o(:)); %# [ sum_r [ sum_c ((o_ij-e_ij)^2/e_ij) ] ] X2 = sum(sum( (o-e).^2 ./ e )); %# p-value needed to reject hNull at the significance level with dof pValue = 1 - chi2cdf(X2, dof); hNull = (pValue > alpha); %# X2 value needed to reject hNull at the significance level with dof %#X2table = chi2inv(1-alpha, dof); %#hNull = (X2table > X2); end

Y un ejemplo para ilustrar:

t = [3 999997 ; 10 999990] [hNull pVal X2] = ChiSquareTest(t, 0.05) hNull = 1 pVal = 0.052203 X2 = 3.7693

Tenga en cuenta que los resultados son diferentes a los suyos porque chisq.test realiza una corrección por defecto, de acuerdo con ?chisq.test

correcto: un indicador lógico que indica si aplicar la corrección de continuidad al calcular el estadístico de prueba para tablas 2x2: una mitad se resta de todos | O - E | diferencias

Alternativamente, si tiene las observaciones reales de los dos eventos en cuestión, puede usar la función CROSSTAB que calcula la tabla de contingencia y devuelve las medidas de Chi2 y p-value:

X = randi([1 2],[1000 1]); Y = randi([1 2],[1000 1]); [t X2 pVal] = crosstab(X,Y) t = 229 247 257 267 X2 = 0.087581 pVal = 0.76728

el equivalente en R sería:

chisq.test(X, Y, correct = FALSE)

Nota: Ambos enfoques (MATLAB) anteriores requieren la Caja de herramientas de estadísticas


Esta función probará la independencia utilizando la estadística de chi-cuadrado de Pearson y la estadística de razón de verosimilitud, junto con el cálculo de los residuos. Sé que esto se puede vectorizar más, pero estoy tratando de mostrar las matemáticas para cada paso.

function independenceTest(data) df = (size(data,1)-1)*(size(data,2)-1); % Mean Degrees of Freedom sd = sqrt(2*df); % Standard Deviation u = nan(size(data)); % Estimated expected frequencies p = nan(size(data)); % Values used to calculate chi-square lr = nan(size(data)); % Values used to calculate likelihood-ratio residuals = nan(size(data)); % Residuals rowTotals = sum(data,1); colTotals = sum(data,2); overallTotal = sum(rowTotals); %% Calculate estimated expected frequencies for r=1:1:size(data,1) for c=1:1:size(data,2) u(r,c) = (rowTotals(c) * colTotals(r)) / overallTotal; end end %% Calculate chi-squared statistic for r=1:1:size(data,1) for c=1:1:size(data,2) p(r,c) = (data(r,c) - u(r,c))^2 / u(r,c); end end chi = sum(sum(p)); % Chi-square statistic %% Calculate likelihood-ratio statistic for r=1:1:size(data,1) for c=1:1:size(data,2) lr(r,c) = data(r,c) * log(data(r,c) / u(r,c)); end end G = 2 * sum(sum(lr)); % Likelihood-Ratio statisitc %% Calculate residuals for r=1:1:size(data,1) for c=1:1:size(data,2) numerator = data(r,c) - u(r,c); denominator = sqrt(u(r,c) * (1 - colTotals(r)/overallTotal) * (1 - rowTotals(c)/overallTotal)); residuals(r,c) = numerator / denominator; end end