0

Matlab で実行する次の実験を念頭に置いており、手順 (3) を実装するための支援を求めています。どんな提案でも大歓迎です。

(1) 確率変数XY両方の一様分布を考える[0,1]

(2)とが独立していると仮定して、 とNの同時分布から実現を導きます (つまり、とは に一様に同時分布することを意味します)。各抽選は で行われます。XYXYXY[0,1]x[0,1][0,1]x[0,1]

(3) ヒルベルト空間充填曲線を使用してドローイン内の各ドローインを変換します[0,1]x[0,1][0,1]ヒルベルト曲線マッピングの下で​​は、ドローインは 内[0,1]x[0,1]の 1 つ (または全射性のために複数) の点の画像になり[0,1]ます。この中から1点選びたいと思います。これを行うMatlabにビルド済みのパッケージはありますか?

ドローのヒルベルト値(曲線の開始点から選択した点までの曲線の長さ)を取得する方法を説明しているので、私が望むことをしているとは思わないこの答えを見つけました

ウィキペディアで、このコードを C 言語 (から(x,y))dで見つけましたが、これも私の質問を満たしていません。

4

4 に答える 4

2

編集 この回答は、ヒルベルト曲線の構築について明示的に尋ねる質問の更新版には対応していません。代わりに、この回答は、全単射写像の構築と一様分布との関係に関する関連する質問に対処します。

あなたの問題は、あまり明確に定義されていません。結果の分布が均一であることだけが必要な場合は、単純に選択することを止めるものは何もありませんf:(X,Y)->XXYが相関しているかどうかに関係なく、結果は一様になります。あなたの投稿から、あなたが実際に望んでいるのは、結果として得られる変換が全単射であるか、機械の精度の制限を考慮して可能な限りそれに近いことであると推測できます。

局所性を維持するのに最適なアルゴリズムが必要でない限り(結果の分布が均一であることは言うまでもなく、全単射であることは明らかに必要ありません)、質問で言及したヒルベルト曲線をわざわざ作成する必要はありません。それらは、他の空間充填曲線と同じように解と関係があり、信じられないほど計算集約的です。

したがって、全単射マッピングを探していると仮定すると、あなたの質問は、[単位] 正方形内の点の集合が [単位] 線分内の点の集合と同じ基数を持つかどうかを尋ねることと同じです。その全単射、つまり 1 対 1 対応を構築する方法。直観は、平方はより高いカーディナリティを持つべきだと言い、カントールはそれを証明しようと3 年を費やし、最終的には正反対のことを証明しました - これらのセットは実際には等数です。彼は自分の発見に非常に驚いたので、次のように書いています。

見えるけど信じられない!

この基準**を満たす、最も一般的に参照される全単射は次のとおりです。andxy10 進形式で表します。つまり、x=0 です。 x1 x2 x3 x4 x5...y=0。 y1 y2 y3 y4 y5...とし、f:(X,Y)->Zz =0 とする。 x1 y1 x2 y2 x3 y3 x4 y4 x5 y5...、つまり、2 つの数値の小数を交互に並べます。全単射の背後にあるアイデアは些細なことですが、厳密な証明にはかなりの事前知識が必要です。

**注意点として、egx = 1/3 = 0.33333...y = 1/5 = 0.199999... = 0.200000...を使用すると、それらに対応する 2 つのシーケンスがあることがわかります:z = 0.313939393939...z = 0.323030303030.... この障害を克服するには、数えられるセットを数えられないセットに追加しても、後者のカーディナリティが変わらないことを証明する必要があります。

実際には、純粋な数学ではなく機械の精度に対処する必要があります。厳密に言えば、両方のセットが実際には有限であり、したがって等数ではないことを意味します (元の数値と同じ精度で結果を格納すると仮定します)。つまり、いくつかの仮定を行うことを余儀なくされ、この場合は と の有効桁数の後半など、いくつかの情報を失うことにxなりyます。つまり、元の変数の精度と比較して、倍精度で結果を格納できる別のデータ型を使用しない限りです。

最後に、Matlab でのサンプル実装:

x = rand();
y = rand();

chars = [num2str(x, '%.17f'); num2str(y, '%.17f')];
z = str2double(['0.' reshape(chars(:,3:end), 1, [])]);

>> cellstr(['x=' num2str(x, '%.17f'); 'y=' num2str(y, '%.17f'); 'z=' num2str(z, '%.17f')])
ans = 
    'x=0.65549803980353738'
    'y=0.10975505072305158'
    'z=0.61505947958500362'
于 2016-05-26T22:29:13.650 に答える
1

You could compute the hilbert curve from f(x,y)=z. Basically it's a hamiltonian path traversal. You can find a good description at Nick's spatial index hilbert curve quadtree blog. Or take a look at monotonic n-ary gray code. I've written an implementation based on Nick's blog in php:http://monstercurves.codeplex.com.

于 2016-05-25T11:42:34.547 に答える
0

私はあなたの最後の点だけに集中します

(3)ヒルベルト空間充填曲線を使用[0,1]x[0,1]して、ドロー イン内の各ドロー インを変換します[0,1]。ヒルベルト曲線マッピングでは、ドロー イン[0,1]x[0,1]は 内の 1 つ (または全射性のために複数) の点のイメージになり[0,1]ます。この中から1点選びたいと思います。これを行うMatlabにビルド済みのパッケージはありますか?

私の知る限り、Matlab にはこれを行うビルド済みのパッケージはありませんが、ウィキペディアのコードは MATLAB から呼び出すことができ、変換ルーチンとゲートウェイ関数を組み合わせるのと同じくらい簡単です。xy2d.cファイル内:

#include "mex.h"

// source: https://en.wikipedia.org/wiki/Hilbert_curve
// rotate/flip a quadrant appropriately
void rot(int n, int *x, int *y, int rx, int ry) {
    if (ry == 0) {
        if (rx == 1) {
            *x = n-1 - *x;
            *y = n-1 - *y;
        }

        //Swap x and y
        int t  = *x;
        *x = *y;
        *y = t;
    }
}

// convert (x,y) to d
int xy2d (int n, int x, int y) {
    int rx, ry, s, d=0;
    for (s=n/2; s>0; s/=2) {
        rx = (x & s) > 0;
        ry = (y & s) > 0;
        d += s * s * ((3 * rx) ^ ry);
        rot(s, &x, &y, rx, ry);
    }
    return d;
}


/* The gateway function */
void mexFunction( int nlhs, mxArray *plhs[],
                  int nrhs, const mxArray *prhs[])
{
    int n;              /* input scalar */
    int x;              /* input scalar */
    int y;              /* input scalar */
    int *d;             /* output scalar */

    /* check for proper number of arguments */
    if(nrhs!=3) {
        mexErrMsgIdAndTxt("MyToolbox:arrayProduct:nrhs","Three inputs required.");
    }
    if(nlhs!=1) {
        mexErrMsgIdAndTxt("MyToolbox:arrayProduct:nlhs","One output required.");
    }

    /* get the value of the scalar inputs  */
    n = mxGetScalar(prhs[0]);
    x = mxGetScalar(prhs[1]);
    y = mxGetScalar(prhs[2]);

    /* create the output */
    plhs[0] = mxCreateDoubleScalar(xy2d(n,x,y));

    /* get a pointer to the output scalar */
    d = mxGetPr(plhs[0]);
}

でコンパイルしmex('xy2d.c')ます。

上記の実装

[...] nセルでnセルに分割された正方形を想定しています。nは 2 の累乗で、左下隅に (0,0)、上隅に ( n -1, n -1) の整数座標があります。右隅。

実際には、マッピングを適用する前に離散化ステップが必要です。すべての離散化問題と同様に、精度を賢く選択することが重要です。以下のスニペットは、すべてをまとめたものです。

close all; clear; clc;

% number of random samples
NSAMPL = 100;

% unit square divided into n-by-n cells
% has to be a power of 2
n = 2^2;

% quantum
d = 1/n;

N = 0:d:1;

% generate random samples
x = rand(1,NSAMPL);
y = rand(1,NSAMPL);

% discretization
bX = floor(x/d);
bY = floor(y/d);

% 2d to 1d mapping
dd = zeros(1,NSAMPL);
for iid = 1:length(dd)
    dd(iid) = xy2d(n, bX(iid), bY(iid));
end


figure;
hold on;
axis equal;

plot(x, y, '.');
plot(repmat([0;1], 1, length(N)), repmat(N, 2, 1), '-r');
plot(repmat(N, 2, 1), repmat([0;1], 1, length(N)), '-r');


figure;
plot(1:NSAMPL, dd);
xlabel('# of sample')
于 2016-06-02T11:03:35.377 に答える