Tích hợp trong Mathematica

Question

Tôi muốn nhận được một giải pháp khác cho một vấn đề mà tôi đã giải quyết "tượng trưng" và thông qua một mô phỏng nhỏ. Bây giờ, tôi muốn biết làm thế nào tôi có thể có sự tích hợp trực tiếp bằng cách sử dụng Mathematica.

Hãy xem xét một mục tiêu được đại diện bởi một đĩa với r = 1, căn giữa tại (0,0) .Tôi muốn mô phỏng xác suất của mình để đạt được mục tiêu ném phi tiêu này.

Bây giờ, tôi không có thành kiến ​​ném chúng, đó là trên trung bình tôi sẽ nhấn trung tâm mu = 0 nhưng sai của tôi là 1.

Xét phối hợp của phi tiêu của tôi vì nó đạt được mục tiêu (hoặc tường :-)) tôi có phân phối sau, 2 Gaussian:

XDistribution : 1/Sqrt[2 \[Pi]\[Sigma]^2] E^(-x^2/(2 \[Sigma]^2)) 

YDistribution : 1/Sqrt[2 \[Pi]\[Sigma]^2] E^(-y^2/(2 \[Sigma]^2)) 

với những 2 phân phối tập trung ở 0 với phương sai bằng = 1, phân phối chung của tôi sẽ trở thành một Gaussian hai biến như:

1/(2 \[Pi]\[Sigma]^2) E^(-((x^2 + y^2)/(2 \[Sigma]^2))) 

Vì vậy, tôi cần biết xác suất của mình để đạt được mục tiêu hoặc xác suất x^2 + y^2 thấp hơn 1.

Tích hợp sau khi chuyển đổi trong hệ tọa độ cực đã cho tôi giải pháp đầu tiên của tôi: .39. Mô phỏng đã xác nhận nó bằng cách sử dụng:

Total@ParallelTable[ 
    If[ 
     EuclideanDistance[{ 
         RandomVariate[NormalDistribution[0, Sqrt[1]]], 
         RandomVariate[NormalDistribution[0, Sqrt[1]]] 
         }, {0, 0}] < 1, 1,0], {1000000}]/1000000 

Tôi cảm thấy có cách thanh lịch hơn để giải quyết vấn đề này bằng cách sử dụng khả năng tích hợp của Mathematica, nhưng không bao giờ phải làm bản đồ ether.

Answer 1

Có một số cách bạn có thể thực hiện việc này.

Cách đơn giản nhất là nên sử dụng NIntegrate như:

JointDistrbution = 1/(2 \[Pi] \[Sigma]^2) E^(-((x^2 + y^2)/(2 \[Sigma]^2))); 
NIntegrate[JointDistrbution /. \[Sigma] -> 1, {y, -1, 1}, 
    {x, -Sqrt[1 - y^2], Sqrt[1 - y^2]}] // Timing 

Out[1]= {0.009625, 0.393469} 

Đây là một cách khác để làm điều đó theo kinh nghiệm (tương tự như ví dụ của bạn ở trên), nhưng chậm hơn rất nhiều so với sử dụng NIntegrate:

(EuclideanDistance[#, {0, 0}] <= 1 & /@ # // Boole // Total)/ 
    Length@# &@RandomVariate[NormalDistribution[0, 1], {10^6, 2}] // 
    N // Timing 

Out[2]= {5.03216, 0.39281} 

Answer 2

Chức năng tích hợp NProbability cũng nhanh:

NProbability[ x^2 + y^2 <= 1, {x, y} \[Distributed] 
BinormalDistribution[{0, 0}, {1, 1}, 0]] // Timing 

hoặc

NProbability[x^2 + y^2 <= 1, x \[Distributed] 
NormalDistribution[0, 1] && y \[Distributed] 
NormalDistribution[0, 1] ] // Timing 

cả hai cung cấp cho {0.031, 0.393469}.

Kể từ khi tổng bình phương của n normals tiêu chuẩn được phân phối ChiSquare[n], bạn sẽ có được một biểu thức hợp lý hơn NProbability[z < 1, z \[Distributed] ChiSquareDistribution[2]] nơi z=x^2+y^2 và x và y được phân phối NormalDistribution[0,1]. Thời gian giống như trên: {0.031, 0.393469}.

EDIT: Thời gian dành cho máy tính 64 bit Core2 Duo T9600 Vista 64 bit với bộ nhớ 8G (MMA 8.0.4). Giải pháp của Yoda trên máy này có thời gian {0.031, 0.393469}.

EDIT 2: Mô phỏng sử dụng ChiSquareDistribution[2] thể được thực hiện như sau:

(data = RandomVariate[ChiSquareDistribution[2], 10^5]; 
    Probability[w <= 1, w \[Distributed] data] // N) // Timing 

mang {0.031, 0.3946}.

EDIT 3: Thông tin chi tiết về timings:

Đối

First@Transpose@Table[Timing@ 
    NProbability[x^2 + y^2 <= 1, {x, y} \[Distributed] 
    BinormalDistribution[{0, 0}, {1, 1}, 0]], {10}] 

tôi nhận được {0.047, 0.031, 0.031, 0.031, 0.031, 0.016, 0.016, 0.031, 0.015, 0.016}

Đối

First@Transpose@Table[Timing@ 
NProbability[x^2 + y^2 <= 1, 
x \[Distributed] NormalDistribution[0, 1] && 
    y \[Distributed] NormalDistribution[0, 1] ], {10}] 

tôi nhận được {0.047, 0.031, 0.032, 0.031, 0.031, 0.016, 0.031, 0.015, 0.016, 0.031}.

Đối

First@Transpose@Table[Timing@ 
NProbability[z < 1, 
z \[Distributed] ChiSquareDistribution[2]], {10}] 

tôi nhận được {0.047, 0.015, 0.016, 0.016, 0.031, 0.015, 0.016, 0.016, 0.015, 0.}.

Đối Yoda của

First@Transpose@Table[Timing@(JointDistrbution = 
    1/(2 \[Pi] \[Sigma]^2) E^(-((x^2 + y^2)/(2 \[Sigma]^2))); 
NIntegrate[ 
    JointDistrbution /. \[Sigma] -> 1, {y, -1, 
    1}, {x, -Sqrt[1 - y^2], Sqrt[1 - y^2]}]), {10}] 

tôi nhận được {0.031, 0.032, 0.015, 0., 0.016, 0., 0.015, 0.016, 0.016, 0.}.

Đối với dự toán thực nghiệm

First@Transpose@Table[Timing@(Probability[w <= 1, 
w \[Distributed] data] // N), {10}] 

tôi đã {0.031, 0.016, 0.016, 0., 0.015, 0.016, 0.015, 0., 0.016, 0.016}.