Area under the curve calculation

[ari]

ID	Time	Concentration
1	0	17
1	0	0
1	0.026	39.5
1	0.031	30.4
1	0.04	29.2
1	0.051	17.8
1	0.092	30.3
1	0.176	71.5
1	0.342	251
1	0.508	291
1	1.015	360
1	2.015	359
1	3.014	617
1	4.01	502
1	5.014	693
1	6.016	672
1	7.019	529
1	8.02	492
1	9.017	590
1	10.013	490
1	11.014	548
1	12.013	443
1	13.013	449
1	14.013	422
1	30.012	295
1	58.019	268
1	120.943	214
1	211.019	185
1	242.995	125
1	272.994	137
1	305.993	156
1	334.998	135
1	366.027	160
1	393.053	127
1	424.996	131
1	456.99	139
1	484.996	120
1	515.981	115
1	544.015	114
1	578.994	133
1	608.986	120
1	638.001	116
1	658.993	122
1	695.088	136
1	726.001	132
1	726.03	237
1	726.051	183
1	726.092	151
1	726.176	155
1	726.342	181
1	727.015	118
1	728.014	72.8
1	729.012	54.6
1	730.021	48.1
2	0	15
2	0	0
2	0.021	27.1
2	0.033	74.9
2	0.052	142
2	0.094	191
2	0.18	197
2	0.346	242
2	0.51	350
2	1.015	442
2	2.015	555
2	3.01	593
2	4.013	540
2	5.016	425
2	6.016	460
2	7.02	610
2	9	514
2	10.001	466
2	11.007	503
2	12.014	415
2	13.024	458
2	14.017	448
2	28.035	44.7
2	56.02	366
2	88.033	278
2	120.026	247
2	215.973	228
2	246.013	186
2	276.978	251
2	306.206	210
2	334.98	190
2	365.073	183
2	393.094	179
2	425.052	217
2	456.056	180
2	484.043	290
2	522.051	317
2	550.059	201
2	582.078	339
2	600.007	203
2	600.045	256
2	600.068	266
2	600.108	183
2	600.191	186
2	600.358	146
2	601.181	81.2
2	602.106	65
2	603.073	31.5
2	604.075	23.2

Hi All,

I need to calculate the area under curve per subject.

I would appreciate any help.

Thanks in advance.

[PGStats]

How do you wish to handle the multiple values at time = 0 (mean, first, last) ?

[ari]

Hi PG,

Thanks for your response.

I would take the last value at time=0.

[PGStats]

Try this:

data want;
do until(last.id);
    set have; by id time;
    if last.time then 
        area = sum(area, mean(concentration, prevConcentration) * 
            range(time, prevTime));
    prevConcentration = concentration;
    prevTime = time;
    end;
keep id area;
run;

Edit: Replaced sum statement with sum function so that the sum doesn't carry from one ID to the next. 

[FreelanceReinh]

Nice implementation. There should be no duplicate time values other than the first for an ID, though.

[ari]

Hi PG,

Thanks for the code.

Is there a way to calculate the AUC in intervals (partial AUC) of 90 days?

Thanks

[PGStats · Accepted solution]

Yes, I think this would do it:

data areas;
startArea = 0;
do until(last.id);
    set have; by id time;
	/* Keep only last observation for any given time value */
    if last.time then do;
		/* Sum a term to trapezoid rule integration */
        area = sum(area, mean(concentration, prevConcentration) * 
            range(time, prevTime));
		/* Look for integration segment limit (multiple of 90) */
        if not missing(prevTime) then 
            if int(time/90) > int(prevTime/90) then do; /* segment limit crossed */
                time90 = 90 * int(time/90);
				/* Interpolate area at segment limit */
                area90 = (prevArea * (time-time90) + area * (time90-prevTime)) / 
                    (time-prevTime) - startArea;
                output;
                startArea = area90;
                end;
        prevConcentration = concentration;
        prevTime = time;
        prevArea = area;
        end;
    end;
keep id time90 area90;
run;

Edit : Added comments. 

[ari]

Thanks PG.

Please explain briefly how the calculation works.

[PGStats]

I added some comments to the code.

[Jack2012]

To be honest, I have not fully understand the mechanism of such algorithm, but it is easier to use and more illustrative compare with PROC IML and post by Rick.

Jack.

[Jack2012]

data have;
infile cards expandtabs truncover;
input ID rho Time 	Concentration;
cards;
1	1 0	17
1	1 0	0
1	1 0.026	39.5
1	1 0.031	30.4
1	1 0.04	29.2
1	1 0.051	17.8
1	1 0.092	30.3
1	1 0.176	71.5
1	1 0.342	251
1	1 0.508	291
1	1 1.015	360
1	1 2.015	359
1	1 3.014	617
1	1 4.01	502
1	1 5.014	693
1	1 6.016	672
1	1 7.019	529
1	1 8.02	492
1	1 9.017	590
1	1 10.013	490
1	1 11.014	548
1	1 12.013	443
1	1 13.013	449
1	2 14.013	422
1	2 30.012	295
1	2 58.019	268
1	2 120.943	214
1	2 211.019	185
1	2 242.995	125
1	2 272.994	137
1	2 305.993	156
1	2 334.998	135
1	2 366.027	160
1	2 393.053	127
1	2 424.996	131
1	2 456.99	139
1	2 484.996	120
1	2 515.981	115
1	2 544.015	114
1	2 578.994	133
1	2 608.986	120
1	2 638.001	116
1	2 658.993	122
1	2 695.088	136
1	2 726.001	132
1	2 726.03	237
1	2 726.051	183
1	2 726.092	151
1	2 726.176	155
1	2 726.342	181
1	2 727.015	118
1	2 728.014	72.8
1	2 729.012	54.6
1	2 730.021	48.1
2	1 0	15
2	1 0	0
2	1 0.021	27.1
2	1 0.033	74.9
2	1 0.052	142
2	1 0.094	191
2	1 0.18	197
2	1 0.346	242
2	1 0.51	350
2	1 1.015	442
2	1 2.015	555
2	1 3.01	593
2	1 4.013	540
2	1 5.016	425
2	1 6.016	460
2	1 7.02	610
2	1 9	514
2	1 10.001	466
2	1 11.007	503
2	1 12.014	415
2	1 13.024	458
2	1 14.017	448
2	1 28.035	44.7
2	1 56.02	366
2	1 88.033	278
2	2 120.026	247
2	2 215.973	228
2	2 246.013	186
2	2 276.978	251
2	2 306.206	210
2	2 334.98	190
2	2 365.073	183
2	2 393.094	179
2	2 425.052	217
2	2 456.056	180
2	2 484.043	290
2	2 522.051	317
2	2 550.059	201
2	2 582.078	339
2	2 600.007	203
2	2 600.045	256
2	2 600.068	266
2	2 600.108	183
2	2 600.191	186
2	2 600.358	146
2	2 601.181	81.2
2	2 602.106	65
2	2 603.073	31.5
2	2 604.075	23.2
;
run;

How should I modify you code to calculate the area by id and rho? Many thanks.

Jack

[Ksharp]

There is trapezoidal rule  which can approximate get the area under curve. Check Rick's blog:

http://blogs.sas.com/content/iml/2011/06/01/the-trapezoidal-rule-of-integration.html

[Jack2012]

Yes, you  are right, I followed this post to calculate the area.

[ari]

Thanks Xia. will check it

[Jack2012]

Dear Keshan,

I have referred to this post by Rick. The question I am facing is how to make use of this by SUBJECT. I mean by use of the PROC IML as used by Rick.

I think there is no BY statement available in PROC IML.

Appreicated for your quick response.

Jack.