Turn on suggestions

Auto-suggest helps you quickly narrow down your search results by suggesting possible matches as you type.

Showing results for

Options

- RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Mute
- Printer Friendly Page

🔒 This topic is **solved** and **locked**.
Need further help from the community? Please
sign in and ask a **new** question.

- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content

Posted 08-16-2020 04:21 PM
(698 views)

Hello

I am running long simulations and using **proc lifereg.** I need to compare the results based on different censoring schemes and different sample size m. One of the problems that I faced is that the Mean squared error (MSE) is not decreasing as sample size increases which is wiered. I thought that maybe I have a convergence problem. So I searched that and the solution that I found on the web is is to fit a log logistic distribution, then use the resulting parameter estimates as initial values in a subsequent fit of a model with the Weibull distribution.

My problem is that I have simulation and in each simulation I get a different data set therefore different parameter estimates using log logistic. How would I solve this problem.

Thank you for your time and help.

1 ACCEPTED SOLUTION

Accepted Solutions

- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content

You didn't copy Ian's code correctly:

`proc lifereg data=AFT outest=est;`

10 REPLIES 10

- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content

Let me answer the MSE question in terms of least squares regression, which is easier to describe and more familiar. The **standard errors** of the parameter estimates should get smaller as the sample size increases. (For linear regression, the standard errors are sqrt(MSE/n).) However, the RMSE = sqrt(MSE) is a measure of the noise component that you add to the linear predictor. It is fixed and does not get smaller with larger sample size.

*> fit a log logistic distribution, then use the resulting parameter estimates as initial values *

OK, I am not familiar with this method, but I understand the main idea of fitting a reduced or associated model.

*> My problem is that I have simulation and in each simulation I get a different data set *

*> therefore different parameter estimates using log logistic. How would I solve this problem?*

I'm not sure what problem you have. Yes, each simulation creates a different sample and therefore different LL param estimates. Use those estimates for the LIFEREG analysis of that sample.

- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content

Thank you Rick

>I'm not sure what problem you have. Yes, each simulation creates a different sample and therefore different LL param estimates. Use those estimates for the LIFEREG analysis of that sample.

The problem is that I can't stop the simulation to copy the results from the first proc lifprog using log-logistic and enter it in the second proc lifprog using weibull. Is there away I can save the results of

proc lifereg data=AFT;

model T1*U_2(0) = Z1/ distribution = llogistic;

run;

Let us call them intercept1 , initial1, scale1 then use them below like this

proc lifereg data=AFT outest=AFT_out covout ;

model T1*U_2(0) = Z1 / itprint dist=WEIBULL

intercept=intercept1 initial= initial1 scale=scale1;

run;

I hope what I am saying makes sense

Thank you.

- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content

I think the trick you are looking for is to use PROC SQL to get the parameter estimates into macro variables. The estimates themselves can be read from the OUTEST data set in LIFEREG.

```
proc lifereg data=AFT outest=est;
model T1*U_2(0) = Z1 / itprint distribution = llogistic;
run;
proc sql noprint;
select intercept format=best20.,
Z1 format=best20.,
_SCALE_ format=best20. into
:p1 TRIMMED, :p2 TRIMMED, :p3 TRIMMED from est;
quit;
%put &=p1 &=p2 &=p3;
```

then in the next model you can specify

`intercept=&p1 initial=&p2 scale=&p3`

as options on the model statement.

- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content

Thank you so much Ian

This is exactly what I wont. However I received an error message which I couldn't figure it out or know how to fix it. Please see the attached file.

Thank you

```
/**AFT wiebul using Progressive TypeII censoring **/
goptions reset=all;
%macro ODSOff(); /* Call prior to BY-group processing */
ods graphics off;
ods exclude all;
ods noresults;
%mend;
%macro ODSOn(); /* Call after BY-group processing */
ods graphics on;
ods exclude none;
ods results;
%mend;
%odsoff;
proc iml;
k=100; seed=0;
B1=0;
scheme=16;n=300; m=280; R=J(m,1,0); do i=1 to m by 14; R[i]=1; end;
******************************;
odds=exp(B1);
lambda=0.5;B0=1;
**** Defining the variables for the Inference part ****;
Beta1=J(k,1,0);
Bias_B1=J(K,1,0);
MSE_B1=J(K,1,0);Var_B1=J(K,1,0);
sq1=J(K,1,0); Lower_bound=J(K,1,0); Upper_bound=J(K,1,0); Test1=J(K,1,0); odd1=J(K,1,0);
MSE_ods=J(K,1,0); Bias_ods=J(K,1,0);
**********************************************************;
Do N1=1 to K; *simulation starts here*;
** Simulation of progressive censoring ****;
hu=J(n,1,0);
Ww=J(m,1,0);X=J(m,1,0);UU1=J(m,1,0);U2=J(m,1,0);T1=J(n,1,0); V=J(n,1,0);
Y1=J(m,1,0);Z1=J(n,1,0);Z2=J(n,1,0);
Ee=J(m,1,0);Vv=J(m,1,0);ll=J(m,1,0);
** Creation of the progressive data from Uniform;
************************************************;
do i=1 to m;
Ww[i]=Uniform(seed);
ll[i]=i;
sum=0;
do j=m-i+1 to m;
sum=(R[j]+sum);
end;
Ee[i]=1/(i+sum);
Vv[i]=Ww[i]**Ee[i];
end;
do i=1 to m;
prod=1;
endloop=m-i+1 ;
do k1=m to endloop by -1;
prod=Vv[k1]*prod;
end;
UU1[i]=1-(prod);** the U's are the required progressively type II
right censored sample from U(0,1);
end;
U_obs=UU1||J(m,1,1);** progressive censored data of size m from U(0,1) with idex=1 to indicate the event occurs;
** I need the data to have a size n, so I need to generate more data from Uniform(0,1) with size of n-m;
UU0=J(n-m,1,0);
do j=1 to n-m;
UU0[j]=Uniform(seed);
end;
U_Cen=UU0||J(n-m,1,0); * These are the censored data with idex 0;
U=J(n,2,0);
U=U_obs//U_cen;
U_1=U[,1];* This is the first column of the matrix U which represents the data from Uniform both(observed + censored);
U_2=U[,2];* This is the index (1:failure OR 0:Censored);
** Simulate SRS data from Weibull **;
do jj=1 to n;
hu[jj]=Uniform(1111);
end;
Z1=quantile('Normal', hu);
C=B0+B1#Z1;
T1=exp(C)#(-log(U_1))##lambda; /*Weibull mixed between progressive and ordinary data**;
/* (1) To use the values of my r.v. in proc lifereg, we need to write them to a SAS data set using CREATE and APPEND
statements. These statements create data set AFT;*/
create AFT var { "U_obs", "U_Cen", "U" , "U_1" ,"U_2","T1", "Z1"};
append;
close AFT;
* (2) Next, the lifereg procedure is called from within IML by using SUBMIT and ENDSUBMIT statements. ;
submit;
proc lifereg data=AFT;
model T1*U_2(0) = Z1/ distribution = llogistic;
run;
proc sql noprint;
select intercept format=best20.,
Z1 format=best20.,
_SCALE_ format=best20. into
:p1 TRIMMED, :p2 TRIMMED, :p3 TRIMMED from est;
quit;
%put &=p1 &=p2 &=p3;
proc lifereg data=AFT outest=AFT_out covout ;
model T1*U_2(0) = Z1 / itprint dist=WEIBULL
intercept=&p1 initial=&p2 scale=&p3;
run;
data new_AFT;
set AFT_out;
keep Z1;
run;
proc transpose data=new_AFT out=idnumber name=Test
prefix=sn;
run;
endsubmit;
* (3) To read the results back to SAS/IML we write "use", "read", and "close";
Use idnumber;
read all var {Sn1 Sn2 Sn3 Sn4};
Close idnumber;
Beta1[N1]=Sn1;
Var_B1[N1]=Sn3;
Use AFT;
read all var {U_2};
Close AFT;
********** Inference ***********;
Beta1[N1]=Sn1;
odd1[N1]=exp(Sn1);
Var_B1[N1]=Sn3;
MSE_B1[N1]=(Beta1[N1]-B1)**2;
MSE_ods[N1]=(odd1[N1]-odds)**2;
end; *end of simulation loop (NN);
%odson
** Estimation of Beta1 ** ;
************************* ;
Beta1_hat=(Beta1[:]);
Bias_B1=abs(Beta1[:]-B1);
Bias_Beta1=round(Bias_B1,0.00001);
MSE_Beta1=round(MSE_B1[:],0.00001);
** Estimation of hazards ration OR Odds ** ;
****************************************** ;
odd1_hat=(odd1[:]);
Biasodd1=abs(odd1_hat-Odds);
Bias_odd1=round(Biasodd1,0.00001);
MSE_Odd1=round(MSE_ods[:],0.00001);
print Bias_Beta1 MSE_Beta1 Bias_odd1 MSE_Odd1;
Quit;
```

- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content

You didn't copy Ian's code correctly:

`proc lifereg data=AFT outest=est;`

- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content

You are right Rick

Thank you so much for your swift reply.

- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content

Thank you Ian

I am so grateful for your trick. It helps me more than you can imagine. AMAZING trick!!

- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content

In addition to the question you asked, you might want to read the article, "Eight tips to make your simulation run faster."

The first tip is to vectorize your IML code by replacing loops with vector operations. For example, in your program, here are two easy loops to replace.

```
/*
sum=0;
do j=m-i+1 to m;
sum=(R[j]+sum);
end;
*/
sum = sum(R[(m-i+1):m]);
...
/*
prod=1;
endloop=m-i+1;
do k1=m to endloop by -1;
prod=Vv[k1]*prod;
end;
*/
prod = prod(Vv[m:(m-i+1)]);
```

- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content

Yes sure sure thank you so much

- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content

You might want to mark Ian's answer as the correct solution.

Registration is open! SAS is returning to Vegas for an AI and analytics experience like no other! Whether you're an executive, manager, end user or SAS partner, SAS Innovate is designed for everyone on your team. Register for just $495 by 12/31/2023.

**If you are interested in speaking, there is still time to submit a session idea. More details are posted on the website. **

Multiple Linear Regression in SAS

Learn how to run multiple linear regression models with and without interactions, presented by SAS user Alex Chaplin.

Find more tutorials on the SAS Users YouTube channel.