In a working paper, Luo (2011) analyzes generating non-normal data, using Fleishman power method. But this approach just can be used for specific skewness and kurtosis. I would like to simulate non-normal data for larger values of skewness and kurtosis (for example, skewness = 3.5 and kurtosis = -4.0). My questions are: 1) how can I do this; 2) are there any R package that can generate such data; 3) how can I generate non-normal data for residuals of the regression model?


Luo, H. (2011). Generation of non-normal data – A study of Fleishman’s power method. [Working paper]. Uppsala University. Retrieved from www.diva-portal.org/smash/get/diva2:407995/FULLTEXT01.pdf

  • $\begingroup$ in addition, actually i wanna generate the non-normal data for residuals of the regression model $\endgroup$ – Nor Hisham Haron May 10 '15 at 0:50
  • $\begingroup$ Please delete your comment, as I embedded it into the question. $\endgroup$ – Aleksandr Blekh May 10 '15 at 1:44
  • $\begingroup$ Kurtosis -4? That's interesting. Are you dealing with complex random numbers? How do you get this? $\endgroup$ – Aksakal May 10 '15 at 1:52
  • $\begingroup$ Please explain which definition of kurtosis you're using there. $\endgroup$ – Glen_b May 10 '15 at 3:00
  • $\begingroup$ As @Aksakal is saying, you might consider another example for skewness and kurtosis. For the most typical definitions (where skewness=kurtosis=0 for a normal distribution), $kurtosis \geq skewness^2 - 2$. $\endgroup$ – Anthony May 11 '15 at 21:23

How to generate non-normal data with specific skewness and kurtosis values?

The Fleishman power method (a.k.a. 3rd order power polynomial) is able to handle a wide range of skewness and kurtosis values, but not all of them. Keep in mind that kurtosis is bounded below by skewness (see Relationship between skew and kurtosis in a sample). The Fleishman method has trouble simulating distributions that are even close to this boundary. A slightly more flexible method is Headrick and Sawilowsky's (1999) 5th order power polynomial method. You can find extensive details about that method and related ones in Headrick's (2009) book. These polynomial methods can work well for many smooth distributions, but they have trouble simulating uniform and most multimodal distributions. For such distributions, you might have to generate some known family (e.g., just generate a uniform with runif()).

Are there any R packages that can generate such data?

Generating the data should be easy using rnorm() and the equations in the cited articles. The hard part is finding the constants to use for specific skewness and kurtosis values. For the Fleishman 3rd order polynomial, Zopluoglu (2011) has R code. For the 5th order power polynomial method, I'm not aware of a reliable solution in R, and I've used Mathematica for that. You can find the Mathematica code here.

How can I generate non-normal data for residuals of the regression model?

Just add your non-normal data to the regression equation. For example, let $E$ be your vector of non-normal residuals:

$Y_i = B_0 + B_1*X_i + E_i$


Headrick, T. C. (2009). Statistical simulation: power method polynomials and other transformations. CRC Press.

Headrick, T. C., & Sawilowsky, S. S. (1999). Simulating correlated multivariate nonnormal distributions: Extending the Fleishman power method. Psychometrika, 64(1), 25-35.

Ruscio, J., & Kaczetow, W. (2008). Simulating multivariate nonnormal data using an iterative algorithm. Multivariate Behavioral Research, 43(3), 355-381.

Zopluoglu, C. (2011). Applications in R: Generating multivariate non-normal variables. http://www.tc.umn.edu/~zoplu001/resim/gennonnormal.pdf.

| cite | improve this answer | |
  • 1
    $\begingroup$ It seems the Zopluoglu paper has moved here: sites.education.miami.edu/zopluoglu/files/2013/12/… $\endgroup$ – jeramy townsley Nov 15 '15 at 6:18
  • 1
    $\begingroup$ I wouldn't take the code written by Zopluoglu too seriously. The approach they use doesn't properly solve the polynomial roots (and no check is made indicate that roots are even found...so impossible combinations are taken as valid). For a more robust and stable approach, try the rValeMaurelli() function from the SimDesign package in R. $\endgroup$ – philchalmers Apr 12 '17 at 18:27

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.