tag:blogger.com,1999:blog-1275149608391671670.post4758907927969063808..comments2019-11-13T01:12:52.087-05:00Comments on SAS and R: Example 7.2: Simulate data from a logistic regressionKen Kleinmanhttp://www.blogger.com/profile/09525118721291529157noreply@blogger.comBlogger24125tag:blogger.com,1999:blog-1275149608391671670.post-65774883344136387532019-08-04T09:26:30.584-04:002019-08-04T09:26:30.584-04:00I'd suggest turning our code into a function t...I'd suggest turning our code into a function then iterating over each of the values in your dataset.Nick Hortonhttps://www.blogger.com/profile/00242216324355342047noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-87261460706993482712019-07-30T12:41:25.147-04:002019-07-30T12:41:25.147-04:00Can someone help me? How do I repeat the simulatio...Can someone help me? How do I repeat the simulation like 30 times having different datasets in one file for example <br />simulation id y x1 x2<br />1 1 1 0 18<br />1 2 0 1 20<br />1 3 0 0 24<br />1 4 1 1 28<br />1 5 0 1 40<br />2 1 1 1 44<br />2 2 0 0 25<br />2 3 1 1 38<br />2 4 1 1 39<br />2 5 1 0 41<br />3 1 1 1 43<br />3 2 1 0 45<br />3 3 1 0 43<br />3 4 0 0 41<br />3 5 0 1 40<br />Seanlovehttps://www.blogger.com/profile/15588123670408061596noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-79710905286488444482019-07-09T14:40:55.407-04:002019-07-09T14:40:55.407-04:00How to simulate a binary response (obese or not ob...How to simulate a binary response (obese or not obese) variable given a distribution of body weight?Anonymousnoreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-51532228599767730722019-04-13T12:27:59.770-04:002019-04-13T12:27:59.770-04:00This is very helpful and informative. I do I intro...This is very helpful and informative. I do I introduce error terms to the ytest?Emmanuelhttps://www.blogger.com/profile/15380787574172947898noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-3657212153716219492018-06-09T21:44:32.702-04:002018-06-09T21:44:32.702-04:00This is very helpful. Thank you for the post!
Is i...This is very helpful. Thank you for the post!<br />Is it possible to modify your program to simulate a case-control data with P(Y=1)=0.5? Many thanks!!Anonymoushttps://www.blogger.com/profile/05770171361554715745noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-24571135764325286882017-11-30T06:05:08.042-05:002017-11-30T06:05:08.042-05:00You can use whatever probabilities you like: the s...You can use whatever probabilities you like: the simulation can be structured to track the scenario of interest.Nick Hortonhttps://www.blogger.com/profile/00242216324355342047noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-78555907759646624732017-11-30T04:23:20.579-05:002017-11-30T04:23:20.579-05:00Thanks, that was what I meant. What I am intereste...Thanks, that was what I meant. What I am interested in is the probabilities of treatment assignment. Is it that I would have one set of probabilities representative of all the datasets or each dataset would have it's own set computed separately? Anonymoushttps://www.blogger.com/profile/17928319076513960795noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-6314509654900545182017-11-27T07:45:40.552-05:002017-11-27T07:45:40.552-05:00Safiya, you could certainly use your dataset as th...Safiya, you could certainly use your dataset as the basis of your simulations and create new Y's using the approach we've described. Is that what you mean by "simulate several logistic regression results"?Nick Hortonhttps://www.blogger.com/profile/00242216324355342047noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-34301181556078238702017-11-27T05:08:11.853-05:002017-11-27T05:08:11.853-05:00Supposing I already have a dataset, can I use same...Supposing I already have a dataset, can I use same to simulate several logistic regression results?Anonymoushttps://www.blogger.com/profile/17928319076513960795noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-55064828816561175322016-10-13T19:36:33.851-04:002016-10-13T19:36:33.851-04:00Mr Ken Kleinman
I used your post how to generate ...Mr Ken Kleinman<br /><br />I used your post how to generate data from a logistic regression in SAS. Its very helpful posts for new users. I use your post of generating data from logistic regression I generate 1000 random numbers, Now I want to replicate this results 100 times, how i can do this. Any suggestions will be pretty helpful. ThanksAnonymousnoreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-11107066336482034012014-06-25T11:46:47.777-04:002014-06-25T11:46:47.777-04:00Possibly. Please contact me privately via e-mail....Possibly. Please contact me privately via e-mail.Ken Kleinmanhttps://www.blogger.com/profile/09525118721291529157noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-59907280684746914802014-06-24T11:12:32.836-04:002014-06-24T11:12:32.836-04:00My primary reason for wanting to do post hoc power...My primary reason for wanting to do post hoc power was to retain the relationships between the variables we have. I'm having trouble simulating variables and accounting for the relationships between them. I have one categorical variable and one continuous one, and I'm interested in looking at the interaction between them. Any suggestions? Anonymoushttps://www.blogger.com/profile/03292349145333273789noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-28780269786967999792014-06-19T11:44:44.235-04:002014-06-19T11:44:44.235-04:00Post-hoc power assessment is fairly controversial ...Post-hoc power assessment is fairly controversial and is frowned upon by many statisticians. But it's not hard to adapt our code to simulate your setting and then to assess power to detect an interaction. I'll write up a new blog post to address this shortly, and thanks for the questionKen Kleinmanhttps://www.blogger.com/profile/09525118721291529157noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-36344703774147871612014-06-18T10:18:19.146-04:002014-06-18T10:18:19.146-04:00I'm interested in doing a post-hoc power test ...I'm interested in doing a post-hoc power test for a logistic regression, but I also have an interaction between continuous and categorical variables. Any suggestions on how to incorporate this into your code? Sorry, I realize this is several years after the original post, just hoping to get some advice as I've only been able to find information on assessing power for an interaction OR for a logistic regression, but not both. Thanks!Anonymoushttps://www.blogger.com/profile/03292349145333273789noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-25893808853843509542012-11-07T13:41:33.607-05:002012-11-07T13:41:33.607-05:00This comment has been removed by the author.Anonymoushttps://www.blogger.com/profile/14530821962623242639noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-28122112617217517162012-11-07T13:40:57.011-05:002012-11-07T13:40:57.011-05:00Hi,
Your program is very helpful, but if I would...Hi, <br /><br />Your program is very helpful, but if I would like to make a restriction of proportion difference like P(Y=1|X=1)-P(Y=1|X=0) = a, is there any way to generate such data controlled by this kind of proportion difference? Any suggestion will be pretty helpful. Thanks!Anonymoushttps://www.blogger.com/profile/14530821962623242639noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-85684944952817720922011-05-25T08:19:49.873-04:002011-05-25T08:19:49.873-04:00With 30 observations, your power is very poor. Th...With 30 observations, your power is very poor. That's the most likely issue. Try simulating and fitting your sample of 30 one thousand times. You should find each parameter with p < .05 more (just slightly more) than 5% of the time. As you increase the sample from 30 to 300, the proportion of rejections should increase-- that's your power increasing.Ken Kleinmanhttps://www.blogger.com/profile/09525118721291529157noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-66969016394289124672011-05-24T22:06:44.365-04:002011-05-24T22:06:44.365-04:00Mr Ken Kleinman
i used your post to generate data ...Mr Ken Kleinman<br />i used your post to generate data logistic. I want to generate x1,x2, and x3. x1 has 2 category (binary), x2 has 3 category, and x3 has normal distribution.<br />When i generate 300, it gave a good model. But when i generate 30, the parameter did't significant. Can you tell me the reason?Anonymousnoreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-89440899119425711552011-03-15T20:47:24.577-04:002011-03-15T20:47:24.577-04:00Right you are. Very reasonable for it to work tha...Right you are. Very reasonable for it to work that way. (I didn't try to make it work that way, obviously, and the documentation doesn't clarify that p can be a vector.) <br /><br />Using your code would be both easier to understand and more "the R way."Ken Kleinmanhttps://www.blogger.com/profile/09525118721291529157noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-30167621669042886492011-03-15T18:16:07.691-04:002011-03-15T18:16:07.691-04:00If I am not mistaken, I think if p is a vector and...If I am not mistaken, I think if p is a vector and rbinom has its first argument is equal to the length of p, this will return a vector with each entry equal to the value in p at that index. <br /><br />rbinom(length(p),1,p)Briannoreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-10523765798495831942011-03-15T09:23:21.351-04:002011-03-15T09:23:21.351-04:00Absolutely. I thought of that last night and wond...Absolutely. I thought of that last night and wondered why didn't write it that way to begin with. I think in R you'd have to use an apply function to allow a different p for each binomial. It would look something like:<br /><br />y = sapply(prob, function(p) {rbinom(1,1,p)})<br /><br />So the binomial part wouls be easier to understand, but the code might be less accessible.Ken Kleinmanhttps://www.blogger.com/profile/09525118721291529157noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-71827143684771556742011-03-14T21:02:41.290-04:002011-03-14T21:02:41.290-04:00I see, that is excellent! The light bulb went on :...I see, that is excellent! The light bulb went on :)<br /><br />I was really thrown by the comparison to the random uniform value. Could you also take a draw from a binomial distribution with p=prob instead?Briannoreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-71631334963424913122011-03-14T20:39:42.782-04:002011-03-14T20:39:42.782-04:00Hi Brian--
ifelse() is a vectorwise logic test. ...Hi Brian--<br /><br />ifelse() is a vectorwise logic test. It tests whether the i element of the runis vector is less than the i element of the prob vector. If it is, then the i element of ytest will be 1, otherwise it will be 0. <br /><br />The prob vector contains the probability of the outcome being 1, given the covariate value and intercept. The ifelse() is like flipping a coin with probability of heads specified by the prob vector.<br /><br />For a sample size calculation via Monte Carlo methods, you would solve prob = exp(intercept)/(1 + exp(intercept)) with prob = known proportion in the population when the covariate=0. You also need to choose the covariate distribution that best mimics your anticipation. Then simulate many times; fit a logistic regression model to each, and see what proportion of times the null is rejected. That's your power. Vary the sample size and beta to find the sample size and beta that achieve the power you need.Ken Kleinmanhttps://www.blogger.com/profile/09525118721291529157noreply@blogger.comtag:blogger.com,1999:blog-1275149608391671670.post-78374789395646011922011-03-14T18:20:14.518-04:002011-03-14T18:20:14.518-04:00Great post, I am confused by ytest = ifelse(runis ...Great post, I am confused by ytest = ifelse(runis < prob,1,0). What is the logic of this step?<br /><br />If you are using this for simulation to determine sample size, would one experiment with the proper intercept (versus zero) to get the probability to average around the expected probability of success in the population?Briannoreply@blogger.com