使用 Schillings 方程计算拟合优度的 P 值答案

【问题标题】：Calculating the P-value for goodness of fits using Schillings equation使用 Schillings 方程计算拟合优度的 P 值
【发布时间】：2015-06-14 17:51:16
【问题描述】：

Nature Methods 最近一篇非常不错的论文给出了一种新方法，用于在不能很好地估计实验误差时计算拟合优度：http://www.nature.com/nmeth/journal/vaop/ncurrent/full/nmeth.3358.html

我正在尝试将其合并到我在 MATLAB 中的拟合例程中，但我的智力低下确实难以正确编码方程式。方程式如下图所示我要测试的是观察到连续补丁（Rn）的概率，因为观察到的最大补丁是 C。

我在 MATLAB 中对这个方程进行了一些微弱的尝试，但结果并不令人印象深刻。我正在努力确定上述方程中的 An，因为当 n>C 和 n更新：我已经更新了代码，因此它不再评估为负概率。但它仍然无法正常工作如果首先创建一个简单的函数来评估 n 的总和

function [ an ] = smalSchill( N, C)
% Small function used to evaluate correlation map statistics

an = 0;
for j = 0:C
    an = an + 2^(N-1-j);
end

end

然后在我的脚本中调用它

C = 10; % observered continoues patch of either positive or negative correlation
N = 100; % length(Q); % Number of observations
g = 1; % dummy variable used to index the P-values
Pv = ones(N,1); % Vector holding the p-values. P = 1 for all n <= C

for n = C+1:N
    An2 = 0;
    for j = 0:C

        if n-1-j <= C
            An2 = An2 + smalSchill(n-1-j,C);
            break
        else
            iter_count = n-1-j-C;
            An2 = An2 + iter_count*smalSchill(C,C); %iter_count            
        end
    end
    Pv2(n) = 1-An2/2^n;
end

P 值开始按预期衰减，但随着 n 的增长，P 值收敛到 1，这是不正确的 - 它们应该收敛到 0。

干杯！

【问题讨论】：

您寻求帮助将论文翻译成 Matlab 代码。为了解决这个问题，需要有人知道或阅读论文，然后与您的代码进行交叉检查。这不是一个真正的编程问题，因此在 SO 上是题外话。
你好唐达。我可能把我的问题表述错了，所以要明确一点：我不要求任何人阅读这篇文章并翻译它。我将文章添加到问题中，给出了我的问题的背景——不是因为你需要阅读它。我正在尝试对出版物中的方程式进行编码，但这样做有困难。任何关于我的代码在方程式方面有什么问题的反馈都将不胜感激 - 我特别想要关于如何避免这么多流控制语句的建议:)
好的，但是如果不深入本文中的细节，我看不出有任何方法可以帮助您解决这个问题。例如。你不解释方程中的符号，给出什么，计算什么等等。
感谢您的反馈 - 我自己也混淆了变量。所以也许我也应该向自己更好地解释一下。下次我一定会更加小心我的问题描述:)

标签： matlab statistics curve-fitting data-fitting

【解决方案1】：

我把一切都搞得太复杂了——并且误解了其中两个变量的含义。难怪我的代码确实有效。但我今天开始工作了

N = 1000 % Length of path investigated
patchL = 1:1:20; % length of continious patch which you are interested in testing.
An = zeros(N+1,length(patchL)); % preallocate memory

for zz = 1:length(patchL)
    x = patchL(zz); 

    for n = 0:N

        if n <= x;
            An(n+1,zz) = 2^n;
        else
            An(n+1,zz) = sum(An(n-x:n,zz));
        end
    end
end

AnDist = An(end,:); % supposing last entry (20) is what you are interested in
Dist = AnDist - [0, AnDist(1:end-1)]; % go from cummulative dist to mdf
Dist = Dist./sum(Dist); % normalizing  the mdf

【讨论】：