3

Trying to fit a chi_square distribution using fitdistr() in R. Documentation on this is here (and not very useful to me): https://stat.ethz.ch/R-manual/R-devel/library/MASS/html/fitdistr.html

Question 1: chi_df below has the following output: 3.85546875 (0.07695236). What is the second number? The Variance or standard deviation?

Question 2: fitdistr generates 'k' defined by the Chi-SQ distribution. How do I fit the data so I get the scaling constant 'A'? I am dumbly using lines 14-17 below. Obviously not good.

Question 3: Is the Chi-SQ distribution only defined for a certain x-range? (Variance is defined as 2K, while mean = k. This must require some constrained x-range... Stats question not programming...)

nnn = 1000;
## Generating a chi-sq distribution
chii <- rchisq(nnn,4, ncp = 0);  
## Plotting Histogram
chi_hist <- hist(chii);   
## Fitting. Gives probability density which must be scaled.
chi_df <- fitdistr(chii,"chi-squared",start=list(df=3)); 
chi_k <- chi_df[[1]][1];

## Plotting a fitted line:
## Spanning x-length of chi-sq data
x_chi_fit <- 1:nnn*((max(chi_hist[[1]][])-min(chi_hist[[1]][]))/nnn);

## Y data using eqn for probability function
y_chi_fit <- (1/(2^(chi_k/2)*gamma(chi_k/2)) * x_chi_fit^(chi_k/2-1) * exp(-x_chi_fit/2));
## Normalizing to the peak of the histogram
y_chi_fit <- y_chi_fit*(max(chi_hist[[2]][]/max(y_chi_fit)));

## Plotting the line
lines(x_chi_fit,y_chi_fit,lwd=2,col="green");

Thanks for your help!

4

1 回答 1

6
  1. 如上所述，?fitdistr说

'"fitdistr"' 类的对象，包含四个组件的列表，... sd：估计的标准误差，

...所以括号中的数字是参数的标准误差。

  1. 尺度参数不需要估计；您需要按直方图箱的宽度进行缩放，或者仅freq=FALSE在绘制直方图时使用。请参阅下面的代码。

  2. 卡方分布是在非负实数上定义的，这是有道理的，因为它是平方标准正态分布（这是一个统计问题，而不是编程问题）。

设置数据：

nnn <- 1000
## ensure reproducibility; not a big deal in this case,
##  but good practice
set.seed(101)
## Generating a chi-sq distribution
chii <- rchisq(nnn,4, ncp = 0)

配件。

library(MASS)
## use method="Brent" based on warning
chi_df <- fitdistr(chii,"chi-squared",start=list(df=3),
                   method="Brent",lower=0.1,upper=100)
chi_k <- chi_df[[1]][1]

（对于它的价值，看起来 print 方法中可能存在一个关于fitdistr何时method="Brent"使用的错误。您也可以使用method="BFGS"并且不需要指定边界......）

直方图

chi_hist <- hist(chii,breaks=50,col="gray")
## scale by N and width of histogram bins
curve(dchisq(x,df=chi_k)*nnn*diff(chi_hist$breaks)[1],
      add=TRUE,col="green")
## or plot histogram already scaled to a density
chi_hist <- hist(chii,breaks=50,col="gray",freq=FALSE)   
curve(dchisq(x,df=chi_k),add=TRUE,col="green")
于 2015-03-04T14:38:34.387 回答