# Data for Question 1
lung=read.table("http://statisticalhorizons.com/wp-content/uploads/LUNG.dat")
names(lung)=c("ID", "INSTIT", "TREAT", "DETECT", "CELL", "STAGE", "TUMOR", "LYMPH", "DISTANT", "OPERATED", "DUR", "SURV")
lung$SURV[lung$SURV==2]=0
lung_extract=lung[c('DUR', 'SURV')]
head(lung_extract, 10)
# a. Construct a survival table similar to the one in the lecture for the
# data. You will show your work within the table – Please use R to obtain the
# survival table (5 Points)
# Load required library
library(survival)
##
## Attaching package: 'survival'
## The following object is masked _by_ '.GlobalEnv':
##
## lung
# Create the Kaplan-Meier survival object
Kaplan.Meier = survfit(Surv(DUR, SURV) ~ 1, data = lung_extract, conf.type = "log-log")
# Get the median survival time
Kaplan.Meier
## Call: survfit(formula = Surv(DUR, SURV) ~ 1, data = lung_extract, conf.type = "log-log")
##
## n events median 0.95LCL 0.95UCL
## [1,] 1032 682 1095 934 1320
# Display the survival table
summary(Kaplan.Meier)
## Call: survfit(formula = Surv(DUR, SURV) ~ 1, data = lung_extract, conf.type = "log-log")
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 0 1032 4 0.996 0.00193 0.990 0.999
## 2 1027 1 0.995 0.00216 0.988 0.998
## 4 1026 1 0.994 0.00237 0.987 0.997
## 6 1025 1 0.993 0.00256 0.986 0.997
## 8 1023 1 0.992 0.00273 0.985 0.996
## 10 1022 1 0.991 0.00290 0.983 0.995
## 12 1021 1 0.990 0.00305 0.982 0.995
## 13 1020 1 0.989 0.00320 0.981 0.994
## 14 1019 2 0.987 0.00348 0.978 0.993
## 16 1017 3 0.984 0.00385 0.975 0.990
## 18 1014 3 0.982 0.00419 0.971 0.988
## 19 1011 1 0.981 0.00430 0.970 0.987
## 21 1010 2 0.979 0.00450 0.968 0.986
## 26 1008 1 0.978 0.00460 0.967 0.985
## 28 1007 1 0.977 0.00470 0.965 0.984
## 29 1006 1 0.976 0.00479 0.964 0.984
## 35 1005 3 0.973 0.00507 0.961 0.981
## 38 1002 3 0.970 0.00532 0.957 0.979
## 39 999 1 0.969 0.00540 0.956 0.978
## 41 998 2 0.967 0.00557 0.954 0.976
## 42 996 2 0.965 0.00572 0.952 0.975
## 43 994 3 0.962 0.00595 0.949 0.972
## 45 991 1 0.961 0.00602 0.947 0.971
## 46 990 2 0.959 0.00616 0.945 0.970
## 49 988 2 0.957 0.00630 0.943 0.968
## 50 986 1 0.956 0.00637 0.942 0.967
## 51 985 1 0.955 0.00643 0.941 0.966
## 60 983 1 0.954 0.00650 0.940 0.966
## 63 982 1 0.953 0.00657 0.939 0.965
## 64 981 1 0.952 0.00663 0.938 0.964
## 65 980 1 0.951 0.00670 0.936 0.963
## 67 979 2 0.950 0.00682 0.934 0.961
## 70 977 2 0.948 0.00694 0.932 0.960
## 72 975 1 0.947 0.00701 0.931 0.959
## 73 974 1 0.946 0.00707 0.930 0.958
## 75 973 1 0.945 0.00712 0.929 0.957
## 76 972 1 0.944 0.00718 0.928 0.956
## 77 971 1 0.943 0.00724 0.927 0.955
## 80 970 1 0.942 0.00730 0.926 0.954
## 86 969 1 0.941 0.00736 0.925 0.954
## 88 968 1 0.940 0.00741 0.923 0.953
## 89 966 1 0.939 0.00747 0.922 0.952
## 94 965 1 0.938 0.00752 0.921 0.951
## 96 963 1 0.937 0.00758 0.920 0.950
## 98 962 1 0.936 0.00763 0.919 0.949
## 102 960 1 0.935 0.00769 0.918 0.948
## 106 959 1 0.934 0.00774 0.917 0.948
## 108 958 1 0.933 0.00779 0.916 0.947
## 110 957 1 0.932 0.00785 0.915 0.946
## 113 955 1 0.931 0.00790 0.914 0.945
## 115 954 3 0.928 0.00805 0.911 0.942
## 117 951 1 0.927 0.00810 0.909 0.941
## 119 950 1 0.926 0.00815 0.908 0.941
## 120 949 2 0.924 0.00825 0.906 0.939
## 122 947 1 0.923 0.00830 0.905 0.938
## 123 946 1 0.922 0.00835 0.904 0.937
## 124 945 1 0.921 0.00840 0.903 0.936
## 125 944 2 0.919 0.00849 0.901 0.934
## 126 942 1 0.918 0.00854 0.900 0.934
## 127 941 3 0.915 0.00868 0.897 0.931
## 128 938 1 0.914 0.00872 0.896 0.930
## 131 937 1 0.913 0.00877 0.895 0.929
## 132 936 1 0.912 0.00881 0.894 0.928
## 134 935 1 0.912 0.00886 0.892 0.927
## 136 934 3 0.909 0.00899 0.889 0.925
## 140 931 1 0.908 0.00903 0.888 0.924
## 141 930 1 0.907 0.00907 0.887 0.923
## 144 929 1 0.906 0.00912 0.886 0.922
## 149 927 1 0.905 0.00916 0.885 0.921
## 150 926 1 0.904 0.00920 0.884 0.920
## 152 925 2 0.902 0.00928 0.882 0.918
## 153 923 3 0.899 0.00941 0.879 0.916
## 154 920 1 0.898 0.00945 0.878 0.915
## 155 919 2 0.896 0.00953 0.876 0.913
## 156 917 1 0.895 0.00957 0.875 0.912
## 161 916 2 0.893 0.00965 0.872 0.910
## 162 914 1 0.892 0.00968 0.871 0.909
## 163 913 2 0.890 0.00976 0.869 0.908
## 165 910 3 0.887 0.00987 0.866 0.905
## 166 907 1 0.886 0.00991 0.865 0.904
## 167 906 1 0.885 0.00995 0.864 0.903
## 169 905 1 0.884 0.00999 0.863 0.902
## 171 904 1 0.883 0.01002 0.862 0.901
## 173 902 2 0.881 0.01010 0.860 0.900
## 175 900 1 0.880 0.01013 0.859 0.899
## 182 899 2 0.878 0.01020 0.857 0.897
## 183 897 3 0.875 0.01031 0.854 0.894
## 184 894 2 0.873 0.01038 0.851 0.892
## 186 892 1 0.872 0.01041 0.850 0.891
## 187 891 1 0.871 0.01045 0.849 0.890
## 190 889 1 0.870 0.01048 0.848 0.890
## 192 888 1 0.869 0.01052 0.847 0.889
## 193 887 1 0.868 0.01055 0.846 0.888
## 194 886 1 0.868 0.01058 0.845 0.887
## 195 885 1 0.867 0.01062 0.844 0.886
## 196 884 1 0.866 0.01065 0.843 0.885
## 203 882 1 0.865 0.01068 0.842 0.884
## 204 881 1 0.864 0.01072 0.841 0.883
## 206 880 1 0.863 0.01075 0.840 0.882
## 207 879 1 0.862 0.01078 0.839 0.881
## 211 878 1 0.861 0.01081 0.838 0.880
## 212 877 2 0.859 0.01088 0.836 0.879
## 213 875 2 0.857 0.01094 0.834 0.877
## 216 871 1 0.856 0.01097 0.833 0.876
## 217 870 1 0.855 0.01100 0.832 0.875
## 220 869 1 0.854 0.01103 0.831 0.874
## 221 868 1 0.853 0.01107 0.830 0.873
## 225 867 1 0.852 0.01110 0.829 0.872
## 226 866 1 0.851 0.01113 0.828 0.871
## 227 864 1 0.850 0.01116 0.826 0.870
## 228 863 1 0.849 0.01119 0.825 0.869
## 230 862 1 0.848 0.01122 0.824 0.868
## 236 859 1 0.847 0.01125 0.823 0.868
## 237 858 1 0.846 0.01128 0.822 0.867
## 238 857 1 0.845 0.01131 0.821 0.866
## 239 856 1 0.844 0.01134 0.820 0.865
## 240 855 1 0.843 0.01137 0.819 0.864
## 241 854 2 0.841 0.01143 0.817 0.862
## 242 852 3 0.838 0.01151 0.814 0.859
## 243 849 1 0.837 0.01154 0.813 0.858
## 245 847 1 0.836 0.01157 0.812 0.857
## 248 846 2 0.834 0.01163 0.810 0.855
## 251 844 1 0.833 0.01166 0.809 0.855
## 253 842 1 0.832 0.01168 0.808 0.854
## 255 841 2 0.830 0.01174 0.806 0.852
## 256 839 1 0.829 0.01177 0.805 0.851
## 259 838 1 0.828 0.01179 0.804 0.850
## 260 837 2 0.826 0.01185 0.801 0.848
## 263 835 1 0.825 0.01188 0.800 0.847
## 264 834 1 0.824 0.01190 0.799 0.846
## 265 833 1 0.823 0.01193 0.798 0.845
## 266 832 1 0.822 0.01196 0.797 0.844
## 274 831 1 0.821 0.01198 0.796 0.843
## 276 830 1 0.820 0.01201 0.795 0.842
## 277 829 1 0.819 0.01204 0.794 0.841
## 279 828 2 0.817 0.01209 0.792 0.840
## 280 826 1 0.816 0.01211 0.791 0.839
## 281 825 1 0.815 0.01214 0.790 0.838
## 283 824 1 0.814 0.01216 0.789 0.837
## 284 823 1 0.813 0.01219 0.788 0.836
## 287 822 1 0.812 0.01222 0.787 0.835
## 294 821 1 0.811 0.01224 0.786 0.834
## 295 819 1 0.810 0.01227 0.785 0.833
## 296 818 1 0.809 0.01229 0.784 0.832
## 299 817 1 0.808 0.01232 0.783 0.831
## 300 815 2 0.806 0.01236 0.781 0.829
## 301 813 1 0.805 0.01239 0.780 0.828
## 303 812 1 0.804 0.01241 0.779 0.827
## 305 811 1 0.803 0.01244 0.778 0.826
## 306 810 1 0.802 0.01246 0.777 0.826
## 307 809 3 0.799 0.01253 0.773 0.823
## 311 806 1 0.798 0.01256 0.772 0.822
## 314 805 1 0.797 0.01258 0.771 0.821
## 315 804 1 0.796 0.01260 0.770 0.820
## 317 803 1 0.795 0.01263 0.769 0.819
## 318 802 1 0.794 0.01265 0.768 0.818
## 319 801 1 0.793 0.01267 0.767 0.817
## 321 799 1 0.792 0.01270 0.766 0.816
## 324 798 1 0.791 0.01272 0.765 0.815
## 326 797 1 0.790 0.01274 0.764 0.814
## 331 796 1 0.789 0.01276 0.763 0.813
## 333 795 1 0.788 0.01279 0.762 0.812
## 334 794 1 0.787 0.01281 0.761 0.811
## 335 793 3 0.784 0.01288 0.758 0.808
## 337 790 1 0.784 0.01290 0.757 0.808
## 338 789 1 0.783 0.01292 0.756 0.807
## 340 788 1 0.782 0.01294 0.755 0.806
## 342 787 1 0.781 0.01296 0.754 0.805
## 347 786 1 0.780 0.01298 0.753 0.804
## 349 785 1 0.779 0.01301 0.752 0.803
## 351 784 1 0.778 0.01303 0.751 0.802
## 356 783 1 0.777 0.01305 0.750 0.801
## 357 782 3 0.774 0.01311 0.747 0.798
## 360 778 1 0.773 0.01313 0.746 0.797
## 361 777 1 0.772 0.01315 0.745 0.796
## 365 776 2 0.770 0.01319 0.742 0.794
## 366 774 1 0.769 0.01321 0.741 0.793
## 368 771 2 0.767 0.01325 0.739 0.791
## 372 768 1 0.766 0.01328 0.738 0.790
## 373 767 1 0.765 0.01330 0.737 0.789
## 374 766 1 0.764 0.01332 0.736 0.789
## 375 765 1 0.763 0.01334 0.735 0.788
## 376 764 1 0.762 0.01336 0.734 0.787
## 379 763 2 0.760 0.01339 0.732 0.785
## 381 761 1 0.759 0.01341 0.731 0.784
## 382 760 1 0.758 0.01343 0.730 0.783
## 384 759 1 0.757 0.01345 0.729 0.782
## 386 758 1 0.756 0.01347 0.728 0.781
## 393 755 3 0.753 0.01353 0.725 0.778
## 395 752 1 0.752 0.01355 0.724 0.777
## 396 750 2 0.750 0.01359 0.722 0.775
## 397 748 1 0.749 0.01361 0.721 0.774
## 399 747 1 0.748 0.01362 0.720 0.773
## 400 746 2 0.746 0.01366 0.718 0.771
## 401 744 1 0.745 0.01368 0.717 0.770
## 403 743 1 0.744 0.01370 0.716 0.769
## 406 740 1 0.743 0.01372 0.715 0.768
## 408 739 1 0.742 0.01373 0.714 0.767
## 417 738 2 0.740 0.01377 0.711 0.765
## 418 736 1 0.739 0.01379 0.710 0.764
## 420 735 1 0.738 0.01381 0.709 0.764
## 422 734 1 0.737 0.01382 0.708 0.763
## 423 733 1 0.736 0.01384 0.707 0.762
## 424 731 1 0.735 0.01386 0.706 0.761
## 428 726 1 0.734 0.01388 0.705 0.760
## 430 725 1 0.733 0.01389 0.704 0.759
## 432 721 1 0.732 0.01391 0.703 0.758
## 433 720 1 0.731 0.01393 0.702 0.757
## 434 719 2 0.728 0.01396 0.700 0.755
## 436 717 1 0.727 0.01398 0.699 0.754
## 438 715 1 0.726 0.01400 0.698 0.753
## 439 713 1 0.725 0.01402 0.697 0.752
## 444 710 1 0.724 0.01403 0.696 0.751
## 446 709 1 0.723 0.01405 0.695 0.750
## 447 707 1 0.722 0.01407 0.694 0.749
## 450 705 1 0.721 0.01409 0.693 0.748
## 451 704 3 0.718 0.01414 0.689 0.745
## 455 701 1 0.717 0.01415 0.688 0.744
## 456 700 2 0.715 0.01419 0.686 0.742
## 457 698 2 0.713 0.01422 0.684 0.740
## 458 696 2 0.711 0.01425 0.682 0.738
## 461 694 1 0.710 0.01427 0.681 0.737
## 464 693 1 0.709 0.01429 0.680 0.736
## 465 692 1 0.708 0.01430 0.679 0.735
## 466 691 1 0.707 0.01432 0.678 0.734
## 470 689 1 0.706 0.01433 0.677 0.733
## 471 688 2 0.704 0.01437 0.675 0.731
## 475 685 1 0.703 0.01438 0.674 0.730
## 477 684 2 0.701 0.01441 0.672 0.728
## 478 681 1 0.700 0.01443 0.670 0.727
## 479 678 1 0.699 0.01444 0.669 0.726
## 480 677 1 0.698 0.01446 0.668 0.725
## 481 676 1 0.697 0.01448 0.667 0.724
## 487 674 1 0.696 0.01449 0.666 0.723
## 488 673 1 0.695 0.01451 0.665 0.722
## 489 672 1 0.694 0.01452 0.664 0.721
## 493 671 1 0.693 0.01454 0.663 0.720
## 495 669 1 0.692 0.01455 0.662 0.719
## 497 668 1 0.690 0.01457 0.661 0.718
## 498 667 1 0.689 0.01458 0.660 0.717
## 500 666 1 0.688 0.01460 0.659 0.716
## 502 665 1 0.687 0.01461 0.658 0.715
## 503 664 1 0.686 0.01463 0.657 0.714
## 504 663 1 0.685 0.01464 0.656 0.713
## 505 662 2 0.683 0.01467 0.654 0.711
## 508 658 1 0.682 0.01468 0.652 0.710
## 509 657 2 0.680 0.01471 0.650 0.708
## 512 655 1 0.679 0.01473 0.649 0.707
## 515 654 3 0.676 0.01477 0.646 0.704
## 517 651 1 0.675 0.01478 0.645 0.703
## 518 650 1 0.674 0.01479 0.644 0.702
## 522 649 2 0.672 0.01482 0.642 0.700
## 531 647 1 0.671 0.01484 0.641 0.699
## 535 646 1 0.670 0.01485 0.640 0.698
## 536 645 1 0.669 0.01486 0.639 0.697
## 542 643 1 0.668 0.01488 0.638 0.696
## 543 642 1 0.667 0.01489 0.637 0.695
## 544 641 3 0.664 0.01493 0.633 0.692
## 545 638 3 0.660 0.01497 0.630 0.689
## 554 635 1 0.659 0.01498 0.629 0.688
## 561 633 1 0.658 0.01499 0.628 0.687
## 563 631 1 0.657 0.01500 0.627 0.686
## 565 630 1 0.656 0.01502 0.626 0.685
## 573 629 1 0.655 0.01503 0.625 0.684
## 575 628 1 0.654 0.01504 0.624 0.683
## 577 626 1 0.653 0.01505 0.623 0.682
## 579 625 1 0.652 0.01506 0.622 0.681
## 582 622 1 0.651 0.01508 0.621 0.680
## 585 621 2 0.649 0.01510 0.618 0.678
## 590 618 1 0.648 0.01511 0.617 0.677
## 592 617 2 0.646 0.01514 0.615 0.675
## 593 615 1 0.645 0.01515 0.614 0.674
## 594 614 1 0.644 0.01516 0.613 0.672
## 596 613 1 0.643 0.01517 0.612 0.671
## 598 612 2 0.641 0.01519 0.610 0.669
## 599 610 1 0.639 0.01521 0.609 0.668
## 601 609 1 0.638 0.01522 0.608 0.667
## 602 608 1 0.637 0.01523 0.607 0.666
## 604 606 1 0.636 0.01524 0.606 0.665
## 608 605 1 0.635 0.01525 0.605 0.664
## 609 604 1 0.634 0.01526 0.603 0.663
## 612 603 1 0.633 0.01527 0.602 0.662
## 617 602 1 0.632 0.01528 0.601 0.661
## 623 601 1 0.631 0.01529 0.600 0.660
## 627 599 1 0.630 0.01530 0.599 0.659
## 631 598 2 0.628 0.01533 0.597 0.657
## 638 596 1 0.627 0.01534 0.596 0.656
## 640 595 1 0.626 0.01535 0.595 0.655
## 641 593 1 0.625 0.01536 0.594 0.654
## 645 592 1 0.624 0.01537 0.593 0.653
## 646 591 1 0.623 0.01538 0.592 0.652
## 648 590 1 0.622 0.01539 0.591 0.651
## 650 589 1 0.620 0.01540 0.590 0.650
## 651 588 2 0.618 0.01542 0.587 0.648
## 652 586 1 0.617 0.01543 0.586 0.647
## 655 585 1 0.616 0.01544 0.585 0.646
## 657 584 1 0.615 0.01545 0.584 0.645
## 662 583 1 0.614 0.01546 0.583 0.644
## 668 581 1 0.613 0.01547 0.582 0.643
## 669 580 1 0.612 0.01547 0.581 0.642
## 677 578 1 0.611 0.01548 0.580 0.641
## 684 577 1 0.610 0.01549 0.579 0.640
## 689 576 1 0.609 0.01550 0.578 0.638
## 690 575 2 0.607 0.01552 0.576 0.636
## 693 573 1 0.606 0.01553 0.575 0.635
## 694 572 3 0.603 0.01556 0.571 0.632
## 696 569 2 0.600 0.01557 0.569 0.630
## 698 567 2 0.598 0.01559 0.567 0.628
## 700 564 1 0.597 0.01560 0.566 0.627
## 717 561 1 0.596 0.01561 0.565 0.626
## 718 560 1 0.595 0.01562 0.564 0.625
## 719 559 1 0.594 0.01562 0.563 0.624
## 738 553 1 0.593 0.01563 0.562 0.623
## 743 549 1 0.592 0.01564 0.561 0.622
## 753 547 1 0.591 0.01565 0.559 0.621
## 755 545 1 0.590 0.01566 0.558 0.620
## 762 544 1 0.589 0.01567 0.557 0.619
## 764 543 1 0.588 0.01568 0.556 0.618
## 766 542 1 0.586 0.01568 0.555 0.616
## 768 541 1 0.585 0.01569 0.554 0.615
## 770 540 1 0.584 0.01570 0.553 0.614
## 771 538 1 0.583 0.01571 0.552 0.613
## 774 536 1 0.582 0.01572 0.551 0.612
## 775 535 1 0.581 0.01573 0.550 0.611
## 778 534 1 0.580 0.01573 0.548 0.610
## 779 533 1 0.579 0.01574 0.547 0.609
## 782 531 1 0.578 0.01575 0.546 0.608
## 793 528 1 0.577 0.01576 0.545 0.607
## 795 527 1 0.576 0.01577 0.544 0.606
## 796 525 1 0.574 0.01577 0.543 0.605
## 797 524 1 0.573 0.01578 0.542 0.604
## 800 523 1 0.572 0.01579 0.541 0.603
## 801 522 1 0.571 0.01580 0.540 0.601
## 804 521 1 0.570 0.01581 0.538 0.600
## 805 520 1 0.569 0.01581 0.537 0.599
## 807 519 1 0.568 0.01582 0.536 0.598
## 811 517 1 0.567 0.01583 0.535 0.597
## 813 514 2 0.565 0.01584 0.533 0.595
## 815 512 2 0.562 0.01586 0.531 0.593
## 816 510 1 0.561 0.01587 0.530 0.592
## 818 509 1 0.560 0.01587 0.528 0.591
## 819 507 1 0.559 0.01588 0.527 0.590
## 824 505 1 0.558 0.01589 0.526 0.588
## 826 504 1 0.557 0.01589 0.525 0.587
## 828 502 1 0.556 0.01590 0.524 0.586
## 836 500 1 0.555 0.01591 0.523 0.585
## 844 498 1 0.554 0.01591 0.522 0.584
## 847 496 2 0.551 0.01593 0.520 0.582
## 848 493 1 0.550 0.01594 0.518 0.581
## 849 492 1 0.549 0.01594 0.517 0.580
## 853 490 1 0.548 0.01595 0.516 0.579
## 856 489 2 0.546 0.01596 0.514 0.576
## 859 487 1 0.545 0.01597 0.513 0.575
## 876 485 1 0.543 0.01598 0.512 0.574
## 882 483 1 0.542 0.01598 0.510 0.573
## 886 481 2 0.540 0.01600 0.508 0.571
## 891 479 1 0.539 0.01600 0.507 0.570
## 895 478 2 0.537 0.01601 0.505 0.568
## 901 476 1 0.536 0.01602 0.504 0.566
## 909 475 1 0.534 0.01603 0.503 0.565
## 921 472 1 0.533 0.01603 0.501 0.564
## 926 471 1 0.532 0.01604 0.500 0.563
## 934 470 2 0.530 0.01605 0.498 0.561
## 951 466 1 0.529 0.01605 0.497 0.560
## 964 464 2 0.526 0.01607 0.495 0.557
## 974 459 2 0.524 0.01608 0.492 0.555
## 977 457 1 0.523 0.01608 0.491 0.554
## 984 455 1 0.522 0.01609 0.490 0.553
## 985 454 1 0.521 0.01609 0.489 0.552
## 991 453 1 0.520 0.01610 0.488 0.551
## 994 452 1 0.518 0.01611 0.486 0.549
## 996 451 1 0.517 0.01611 0.485 0.548
## 999 450 2 0.515 0.01612 0.483 0.546
## 1038 441 1 0.514 0.01613 0.482 0.545
## 1043 440 1 0.513 0.01613 0.481 0.544
## 1049 438 1 0.512 0.01614 0.479 0.543
## 1051 437 1 0.510 0.01614 0.478 0.541
## 1052 436 1 0.509 0.01615 0.477 0.540
## 1066 434 1 0.508 0.01615 0.476 0.539
## 1077 433 1 0.507 0.01616 0.475 0.538
## 1085 432 2 0.504 0.01617 0.472 0.536
## 1091 429 1 0.503 0.01617 0.471 0.535
## 1093 428 1 0.502 0.01618 0.470 0.533
## 1095 427 2 0.500 0.01619 0.468 0.531
## 1104 424 1 0.499 0.01619 0.466 0.530
## 1105 423 1 0.497 0.01620 0.465 0.529
## 1114 422 1 0.496 0.01620 0.464 0.528
## 1116 421 1 0.495 0.01621 0.463 0.526
## 1118 420 1 0.494 0.01621 0.462 0.525
## 1119 419 1 0.493 0.01622 0.461 0.524
## 1143 415 1 0.492 0.01622 0.459 0.523
## 1159 413 1 0.490 0.01622 0.458 0.522
## 1161 412 1 0.489 0.01623 0.457 0.520
## 1162 411 1 0.488 0.01623 0.456 0.519
## 1172 409 1 0.487 0.01624 0.455 0.518
## 1178 408 1 0.486 0.01624 0.453 0.517
## 1180 407 1 0.484 0.01624 0.452 0.516
## 1191 406 1 0.483 0.01625 0.451 0.515
## 1198 404 1 0.482 0.01625 0.450 0.513
## 1208 402 1 0.481 0.01626 0.449 0.512
## 1211 401 1 0.480 0.01626 0.447 0.511
## 1216 400 1 0.478 0.01626 0.446 0.510
## 1218 399 1 0.477 0.01627 0.445 0.509
## 1220 397 1 0.476 0.01627 0.444 0.507
## 1234 395 1 0.475 0.01627 0.443 0.506
## 1284 389 1 0.474 0.01628 0.441 0.505
## 1290 388 1 0.472 0.01628 0.440 0.504
## 1292 387 1 0.471 0.01628 0.439 0.503
## 1295 386 1 0.470 0.01629 0.438 0.501
## 1311 382 1 0.469 0.01629 0.436 0.500
## 1320 381 1 0.467 0.01629 0.435 0.499
## 1323 380 1 0.466 0.01630 0.434 0.498
## 1325 379 1 0.465 0.01630 0.433 0.497
## 1329 377 1 0.464 0.01630 0.431 0.495
## 1330 376 1 0.462 0.01631 0.430 0.494
## 1331 375 1 0.461 0.01631 0.429 0.493
## 1336 373 1 0.460 0.01631 0.428 0.492
## 1356 366 1 0.459 0.01632 0.427 0.490
## 1359 364 1 0.458 0.01632 0.425 0.489
## 1364 363 1 0.456 0.01633 0.424 0.488
## 1370 362 1 0.455 0.01633 0.423 0.487
## 1371 361 1 0.454 0.01633 0.421 0.485
## 1391 359 1 0.452 0.01634 0.420 0.484
## 1402 358 1 0.451 0.01634 0.419 0.483
## 1429 357 1 0.450 0.01634 0.418 0.482
## 1431 356 1 0.449 0.01634 0.416 0.480
## 1441 355 1 0.447 0.01635 0.415 0.479
## 1457 351 1 0.446 0.01635 0.414 0.478
## 1462 350 1 0.445 0.01635 0.413 0.477
## 1488 346 1 0.444 0.01636 0.411 0.475
## 1490 345 1 0.442 0.01636 0.410 0.474
## 1493 343 1 0.441 0.01636 0.409 0.473
## 1494 342 1 0.440 0.01637 0.407 0.471
## 1500 341 1 0.438 0.01637 0.406 0.470
## 1505 340 1 0.437 0.01637 0.405 0.469
## 1511 339 1 0.436 0.01637 0.404 0.468
## 1515 338 1 0.435 0.01638 0.402 0.466
## 1535 336 1 0.433 0.01638 0.401 0.465
## 1542 335 1 0.432 0.01638 0.400 0.464
## 1548 333 1 0.431 0.01638 0.398 0.463
## 1550 332 1 0.429 0.01638 0.397 0.461
## 1559 330 1 0.428 0.01639 0.396 0.460
## 1569 329 1 0.427 0.01639 0.394 0.459
## 1578 327 1 0.425 0.01639 0.393 0.457
## 1583 325 1 0.424 0.01639 0.392 0.456
## 1600 320 1 0.423 0.01639 0.391 0.455
## 1606 318 1 0.421 0.01640 0.389 0.453
## 1608 317 1 0.420 0.01640 0.388 0.452
## 1616 316 1 0.419 0.01640 0.387 0.451
## 1617 315 1 0.418 0.01640 0.385 0.449
## 1627 314 1 0.416 0.01640 0.384 0.448
## 1630 313 1 0.415 0.01640 0.383 0.447
## 1638 311 1 0.414 0.01641 0.381 0.445
## 1640 310 1 0.412 0.01641 0.380 0.444
## 1643 309 1 0.411 0.01641 0.379 0.443
## 1644 308 1 0.410 0.01641 0.377 0.441
## 1663 305 1 0.408 0.01641 0.376 0.440
## 1665 304 1 0.407 0.01641 0.375 0.439
## 1668 303 1 0.405 0.01641 0.373 0.437
## 1682 301 1 0.404 0.01641 0.372 0.436
## 1708 295 1 0.403 0.01641 0.371 0.435
## 1713 294 1 0.401 0.01641 0.369 0.433
## 1717 293 1 0.400 0.01642 0.368 0.432
## 1720 291 1 0.399 0.01642 0.366 0.431
## 1731 287 1 0.397 0.01642 0.365 0.429
## 1735 285 1 0.396 0.01642 0.364 0.428
## 1760 282 2 0.393 0.01642 0.361 0.425
## 1769 278 1 0.392 0.01642 0.359 0.424
## 1781 274 1 0.390 0.01643 0.358 0.422
## 1786 272 1 0.389 0.01643 0.357 0.421
## 1787 270 1 0.387 0.01643 0.355 0.419
## 1799 268 1 0.386 0.01643 0.354 0.418
## 1806 267 1 0.384 0.01644 0.352 0.417
## 1808 263 1 0.383 0.01644 0.351 0.415
## 1815 261 1 0.382 0.01644 0.349 0.414
## 1819 260 1 0.380 0.01644 0.348 0.412
## 1831 257 1 0.379 0.01644 0.346 0.411
## 1836 253 1 0.377 0.01645 0.345 0.409
## 1844 250 1 0.376 0.01645 0.343 0.408
## 1845 248 1 0.374 0.01645 0.342 0.406
## 1849 246 1 0.373 0.01646 0.340 0.405
## 1850 245 1 0.371 0.01646 0.339 0.403
## 1851 244 1 0.369 0.01646 0.337 0.402
## 1856 240 1 0.368 0.01647 0.336 0.400
## 1862 239 1 0.366 0.01647 0.334 0.399
## 1871 237 1 0.365 0.01647 0.333 0.397
## 1875 234 1 0.363 0.01647 0.331 0.396
## 1881 233 1 0.362 0.01648 0.330 0.394
## 1886 230 1 0.360 0.01648 0.328 0.392
## 1890 229 1 0.359 0.01648 0.326 0.391
## 1898 228 1 0.357 0.01649 0.325 0.389
## 1909 227 1 0.355 0.01649 0.323 0.388
## 1921 224 1 0.354 0.01649 0.322 0.386
## 1922 223 1 0.352 0.01649 0.320 0.385
## 1936 220 1 0.351 0.01650 0.318 0.383
## 1937 219 1 0.349 0.01650 0.317 0.381
## 1939 218 1 0.347 0.01650 0.315 0.380
## 1943 217 1 0.346 0.01650 0.314 0.378
## 1944 215 1 0.344 0.01650 0.312 0.377
## 1945 214 1 0.343 0.01650 0.310 0.375
## 1951 211 1 0.341 0.01651 0.309 0.373
## 1960 209 1 0.339 0.01651 0.307 0.372
## 1963 208 1 0.338 0.01651 0.306 0.370
## 1965 207 1 0.336 0.01651 0.304 0.369
## 1980 206 1 0.335 0.01651 0.302 0.367
## 1985 205 1 0.333 0.01651 0.301 0.365
## 1987 204 1 0.331 0.01651 0.299 0.364
## 1992 203 1 0.330 0.01651 0.297 0.362
## 2001 200 1 0.328 0.01651 0.296 0.360
## 2002 199 1 0.326 0.01651 0.294 0.359
## 2006 198 1 0.325 0.01651 0.293 0.357
## 2017 194 1 0.323 0.01650 0.291 0.355
## 2028 192 1 0.321 0.01650 0.289 0.354
## 2041 189 1 0.320 0.01650 0.288 0.352
## 2042 188 1 0.318 0.01650 0.286 0.350
## 2043 187 1 0.316 0.01650 0.284 0.349
## 2044 186 1 0.315 0.01650 0.282 0.347
## 2057 182 1 0.313 0.01650 0.281 0.345
## 2069 179 1 0.311 0.01650 0.279 0.344
## 2073 178 1 0.309 0.01650 0.277 0.342
## 2074 177 1 0.308 0.01650 0.276 0.340
## 2084 176 1 0.306 0.01650 0.274 0.338
## 2107 173 1 0.304 0.01650 0.272 0.337
## 2110 172 1 0.302 0.01650 0.270 0.335
## 2125 170 1 0.300 0.01649 0.269 0.333
## 2148 164 1 0.299 0.01650 0.267 0.331
## 2150 163 1 0.297 0.01650 0.265 0.329
## 2155 161 1 0.295 0.01650 0.263 0.328
## 2168 158 1 0.293 0.01650 0.261 0.326
## 2169 157 1 0.291 0.01650 0.259 0.324
## 2171 156 1 0.289 0.01650 0.257 0.322
## 2174 155 1 0.288 0.01650 0.256 0.320
## 2178 153 1 0.286 0.01649 0.254 0.318
## 2187 151 1 0.284 0.01649 0.252 0.316
## 2189 150 1 0.282 0.01649 0.250 0.315
## 2190 149 2 0.278 0.01649 0.246 0.311
## 2191 146 1 0.276 0.01648 0.244 0.309
## 2200 143 1 0.274 0.01648 0.242 0.307
## 2210 141 1 0.272 0.01648 0.240 0.305
## 2222 136 1 0.270 0.01648 0.238 0.303
## 2270 133 1 0.268 0.01648 0.236 0.301
## 2287 129 1 0.266 0.01648 0.234 0.299
## 2288 128 1 0.264 0.01648 0.232 0.297
## 2301 124 1 0.262 0.01649 0.230 0.295
## 2311 121 1 0.260 0.01649 0.228 0.293
## 2341 117 1 0.258 0.01650 0.226 0.290
## 2353 115 1 0.255 0.01651 0.224 0.288
## 2356 114 1 0.253 0.01651 0.221 0.286
## 2363 111 1 0.251 0.01652 0.219 0.284
## 2373 110 1 0.249 0.01653 0.217 0.281
## 2376 109 1 0.246 0.01653 0.215 0.279
## 2382 108 1 0.244 0.01654 0.212 0.277
## 2398 106 1 0.242 0.01654 0.210 0.275
## 2407 105 1 0.239 0.01654 0.208 0.272
## 2414 104 1 0.237 0.01654 0.205 0.270
## 2419 102 1 0.235 0.01654 0.203 0.268
## 2450 97 1 0.232 0.01655 0.201 0.265
## 2468 93 1 0.230 0.01656 0.198 0.263
## 2480 92 1 0.227 0.01657 0.196 0.260
## 2495 89 1 0.225 0.01657 0.193 0.258
## 2500 88 1 0.222 0.01658 0.191 0.255
## 2521 83 1 0.220 0.01660 0.188 0.253
## 2523 82 1 0.217 0.01661 0.185 0.250
## 2532 80 1 0.214 0.01662 0.183 0.248
## 2536 79 1 0.211 0.01663 0.180 0.245
## 2539 78 1 0.209 0.01664 0.177 0.242
## 2555 75 1 0.206 0.01665 0.174 0.239
## 2560 74 1 0.203 0.01665 0.172 0.237
## 2691 60 1 0.200 0.01672 0.168 0.233
## 2722 56 1 0.196 0.01679 0.164 0.230
## 2772 52 1 0.192 0.01689 0.161 0.227
## 2823 47 1 0.188 0.01702 0.156 0.223
## 2955 37 1 0.183 0.01730 0.151 0.218
## 3053 30 1 0.177 0.01777 0.144 0.213
## 3156 25 1 0.170 0.01842 0.136 0.208
# Kaplan-Meier Survival Model
# Model:
# The analysis uses the Kaplan-Meier non-parametric estimator to analyze
# survival data:
# survfit(Surv(DUR, SURV) ~ 1, data = lung_extract, conf.type = "log-log")
# Model Elements:
# Outcome variable: Survival time (DUR) and event indicator (SURV)
# Estimator: Product-limit (Kaplan-Meier) estimator
# Confidence interval method: Log-log transformation
# Assumptions:
# Independent observations: Each patient's survival time is independent of
# other patients
# Non-informative censoring: The censoring process is independent of the
# survival process
# Constant hazard between observed event times: The risk of event remains
# constant between observed times
# Equivalence of subjects: All subjects at risk at each time point have the
# same risk of the event.
# Metric Value
# Total observations 1032
# Events (deaths) 682
# Median survival 1095 days
# 95% CI for median 934 to 1320 days
# Conclusion:
# Based on the Kaplan-Meier analysis:
# Median survival time is 1,095 days, with 95% CI from 934 to 1320 days.
# The 1-year survival probability is approximately 0.81 (estimated from time
# points near 365 days)
# The 2-year survival probability is approximately 0.65 (estimated from time
# points near 730 days)
# The survival probability decreases over time, as expected in a population
# with lung cancer
# The precision of estimates (confidence intervals)
# b. Plot a survival curve for the probabilities you generate in part a (5 Points)
# Load required library
library(survival)
# Create Kaplan-Meier survival object
Kaplan.Meier = survfit(Surv(DUR, SURV) ~ 1, data = lung_extract, conf.type = "log-log")
# Plot the survival curve
plot(Kaplan.Meier,
main = "Survival Curve - Kaplan Meier",
xaxt = "n",
xlab = "time(days)",
ylab = "Cumulative Survival")
# Add custom x-axis with appropriate intervals
axis(1, at = seq(0, 4000, by = 500), las = 2)
# Add grid for better readability
grid()
# Add confidence intervals if desired
lines(Kaplan.Meier$time, Kaplan.Meier$upper, type = "s", col = "red", lty = 2)
lines(Kaplan.Meier$time, Kaplan.Meier$lower, type = "s", col = "red", lty = 2)
# Add a legend
legend("topright",
legend = c("Survival curve", "95% Confidence Interval"),
col = c("black", "red"),
lty = c(1, 2),
lwd = c(1, 1))
# The survival curve shows:
# At time 0, survival probability is 100% (everyone is alive at the start)
# There's a steady decline in survival probability over the first 2,500 days
# By around 3,000 days, the survival probability plateaus at approximately 20%
# After the 3,000-day mark, the curve flattens, suggesting that patients who
# survive beyond this point have a relatively stable outlook
# The widening confidence interval toward the end indicates increased
# uncertainty in the estimates due to fewer remaining subjects
# This curve is specifically analyzing lung cancer mortality, where status=1
# indicates death from the disease and status=0 indicates the patient was
# either still alive or died from other causes (both treated as censoring in
# this analysis).
# The plot effectively communicates the long-term survival outlook for
# patients in this lung dataset, showing that approximately 20% of patients
# survive beyond 8-10 years (3,000 plus days).
2.You are studying three different new drugs that may help slow the progress of La Traviata disease which compels people to sing opera until they exhaust themselves and die. Do the following:
Data for part 2 of this exercise in the attached R script. Data format is Group#, time event (0=no event, 1=event)
# Data for Question 2
Traviata=data.frame(
group=c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3),
time=c(681, 602, 996, 1162, 833, 477, 630, 596, 226, 699, 811, 530, 482, 367, 118, 83, 76, 104, 109, 72, 87, 162, 94, 30,
26, 22, 49, 74, 122, 86, 66, 92, 109, 255, 1, 107, 110, 232, 2569, 2506, 2409, 2218, 1857, 1829, 1562, 1470, 1363,
1030, 1860, 1258, 2246, 1870, 1799, 1709, 1674, 1568, 1527, 1324, 1957, 1932, 1847, 1848, 1850, 1843, 1535, 1447,
1384, 914, 2204, 1063, 481, 605, 641, 390, 288, 421, 1379, 1748, 486, 448, 272, 1074, 1381, 1410, 1353, 1480, 435,
248, 1704, 1411, 219, 606, 2640, 2430, 2252, 2140, 2133, 1738, 2631, 2524, 1845, 1936, 1845, 422, 162, 84, 100, 212,
47, 242, 456, 268, 318, 732, 467, 947, 390, 183, 105, 115, 164, 693, 120, 80, 677, 64, 168, 874, 616, 157, 625, 48,
273, 163, 376, 113, 363),
event=c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1)
)
# Convert group to factor and label it appropriately
Traviata$group <- factor(Traviata$group, levels = c(1, 2, 3), labels = c("Drug 1", "Drug 2", "Drug 3"))
# Print summary statistics by group
cat("Summary statistics by drug group:\n")
## Summary statistics by drug group:
summary_stats <- table(Traviata$group, Traviata$event)
colnames(summary_stats) <- c("Censored", "Events")
print(summary_stats)
##
## Censored Events
## Drug 1 14 24
## Drug 2 29 25
## Drug 3 11 34
# Calculate total counts and event rates
group_counts <- table(Traviata$group)
event_rates <- summary_stats[, "Events"] / group_counts
cat("\nEvent rates by group:\n")
##
## Event rates by group:
print(data.frame(Group = rownames(summary_stats), Count = group_counts,
Events = summary_stats[, "Events"],
EventRate = event_rates))
## Group Count.Var1 Count.Freq Events EventRate.Var1 EventRate.Freq
## Drug 1 Drug 1 Drug 1 38 24 Drug 1 0.6315789
## Drug 2 Drug 2 Drug 2 54 25 Drug 2 0.4629630
## Drug 3 Drug 3 Drug 3 45 34 Drug 3 0.7555556
# MODEL AND ASSUMPTIONS
# Kaplan-Meier model for survival analysis with the following assumptions
# Censoring is non-informative, censored patients have same survival
# prospects as those who continue
# Survival probabilities depend only on time since entry into study
# Patients recruited early and late in the study have similar survival
# prospects
# Events happen at the specified times
# Independence between observations
# Create survival object
surv_obj <- Surv(time = Traviata$time, event = Traviata$event)
# Fit Kaplan-Meier curves for each group
km_fit <- survfit(surv_obj ~ group, data = Traviata)
print(km_fit)
## Call: survfit(formula = surv_obj ~ group, data = Traviata)
##
## n events median 0.95LCL 0.95UCL
## group=Drug 1 38 24 114 94 NA
## group=Drug 2 54 25 2204 1411 NA
## group=Drug 3 45 34 376 212 693
# 3. MODEL ELEMENTS ANNOTATION
# Survival object (Surv): Time-to-event data with censoring indicator
# Kaplan-Meier estimator (survfit): Non-parametric method to estimate
# survival function
# Log-rank test (survdiff): Non-parametric test to compare survival distributions
# Pairwise comparisons: Specific tests between pairs of drugs
# Group factor: Treatment group (Drug 1, Drug 2, Drug 3)
# Time: Time until event or censoring (in days)
# Event: Binary indicator (1 = event occurred, 0 = censored)
# A. SURVIVAL PLOT FOR ALL THREE DRUGS
cat("\n--- 4A. SURVIVAL PLOT FOR ALL THREE DRUGS ---\n")
##
## --- 4A. SURVIVAL PLOT FOR ALL THREE DRUGS ---
# Create survival plot
plot(km_fit,
main = "Kaplan-Meier Survival Curves for La Traviata Disease",
xlab = "Time (days)",
ylab = "Survival Probability",
col = c("blue", "red", "green"),
lwd = 2)
# Add legend
legend("topright",
legend = c("Drug 1", "Drug 2", "Drug 3"),
col = c("blue", "red", "green"),
lwd = 2,
title = "Treatment Group")
# Drug 2 (red line) shows the best survival outcomes:
# Maintains the highest survival probability throughout the observation period
# Has approximately 42% survival probability at the end of the study (around 2,500 days)
# Shows a more gradual decline in survival compared to the other treatments
# Drug 3 (green line) shows intermediate effectiveness:
# Initially performs better than Drug 1
# Ends with approximately 22% survival probability
# Shows a steep initial decline within the first 500 days
# Drug 1 (blue line) appears least effective:
# Shows the most rapid decline in survival probability within the first 300 days
# Plateaus at about 35% survival probability around 500 days
# Data collection for this group appears to end earlier (around 1,000 days)
# General timeline patterns:
# The first 500 days show the steepest decline for all treatment groups
# After about 1,000 days, the curves begin to plateau, suggesting lower
# mortality risk for patients who survive beyond this point
# This plot suggests Drug 2 would be the most promising treatment
# option based on long-term survival outcomes for patients with La Traviata Disease.
#b. Test to see if overall there is an effect of any of the drugs on survival taken as a global set (25 Points).
#TEST FOR OVERALL EFFECT OF DRUGS ON SURVIVAL
cat("\n--- 4B. TEST FOR OVERALL EFFECT OF DRUGS ON SURVIVAL ---\n")
##
## --- 4B. TEST FOR OVERALL EFFECT OF DRUGS ON SURVIVAL ---
# Null hypothesis (H0): There is no difference in survival between any of the drug groups
# Alternative hypothesis (H1): At least one drug has a different survival distribution
# Perform log-rank test
log_rank_test <- survdiff(surv_obj ~ group, data = Traviata)
cat("Log-rank test for differences in survival across all three drugs:\n")
## Log-rank test for differences in survival across all three drugs:
print(log_rank_test)
## Call:
## survdiff(formula = surv_obj ~ group, data = Traviata)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## group=Drug 1 38 24 12.3 11.07 13.66
## group=Drug 2 54 25 46.0 9.58 22.50
## group=Drug 3 45 34 24.7 3.51 5.04
##
## Chisq= 25.7 on 2 degrees of freedom, p= 3e-06
# Calculate p-value
p_value <- 1 - pchisq(log_rank_test$chisq, df = length(levels(Traviata$group)) - 1)
cat("\nP-value:", format.pval(p_value, digits = 6), "\n")
##
## P-value: 2.66883e-06
# Interpretation
alpha <- 0.05
cat("\nHypothesis Test:\n")
##
## Hypothesis Test:
cat("H₀: Survival curves are equal across all 3 groups.\n")
## H₀: Survival curves are equal across all 3 groups.
cat("H₁: At least one survival curve is different.\n")
## H₁: At least one survival curve is different.
cat("Test statistic: χ² =", round(log_rank_test$chisq, 2), "with df =", length(levels(Traviata$group)) - 1, "\n")
## Test statistic: χ² = 25.67 with df = 2
cat("Decision rule: Reject H₀ if p-value <", alpha, "\n")
## Decision rule: Reject H₀ if p-value < 0.05
if(p_value < alpha) {
cat("Result: Reject H₀ at alpha =", alpha, "\n")
cat("Conclusion: There is a statistically significant difference in survival between at least two drug groups.\n")
} else {
cat("Result: Fail to reject H₀ at alpha =", alpha, "\n")
cat("Conclusion: There is no statistically significant difference in survival between the drug groups.\n")
}
## Result: Reject H₀ at alpha = 0.05
## Conclusion: There is a statistically significant difference in survival between at least two drug groups.
# c. Compare the survival curves for each of the three drugs with each other
# (three comparisons) and see if any if the curves are different from each
# other (10 Points).
# PAIRWISE COMPARISONS BETWEEN DRUGS
cat("\n--- 4C. PAIRWISE COMPARISONS BETWEEN DRUGS ---\n")
##
## --- 4C. PAIRWISE COMPARISONS BETWEEN DRUGS ---
# Create subsets for pairwise comparisons
Traviata_1vs2 <- Traviata[Traviata$group %in% c("Drug 1", "Drug 2"), ]
Traviata_1vs3 <- Traviata[Traviata$group %in% c("Drug 1", "Drug 3"), ]
Traviata_2vs3 <- Traviata[Traviata$group %in% c("Drug 2", "Drug 3"), ]
# Set significance level and apply Bonferroni correction for multiple tests
alpha <- 0.05
alpha_corrected <- alpha / 3 # Bonferroni correction for 3 comparisons
cat("Using Bonferroni-corrected significance level:", alpha_corrected, "\n\n")
## Using Bonferroni-corrected significance level: 0.01666667
# 1. Drug 1 vs Drug 2
cat("--- Drug 1 vs Drug 2 ---\n")
## --- Drug 1 vs Drug 2 ---
# Create new survival object for this subset
surv_obj_1vs2 <- Surv(time = Traviata_1vs2$time, event = Traviata_1vs2$event)
# Perform log-rank test
drug1_vs_drug2 <- survdiff(surv_obj_1vs2 ~ group, data = Traviata_1vs2)
print(drug1_vs_drug2)
## Call:
## survdiff(formula = surv_obj_1vs2 ~ group, data = Traviata_1vs2)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## group=Drug 1 38 24 10.3 18.29 26.1
## group=Drug 2 54 25 38.7 4.86 26.1
##
## Chisq= 26.1 on 1 degrees of freedom, p= 3e-07
p_value_1vs2 <- 1 - pchisq(drug1_vs_drug2$chisq, df = 1)
cat("P-value:", format.pval(p_value_1vs2, digits = 6), "\n")
## P-value: 3.27492e-07
cat("Significant at α =", alpha_corrected, ":", p_value_1vs2 < alpha_corrected, "\n\n")
## Significant at α = 0.01666667 : TRUE
# 2. Drug 1 vs Drug 3
cat("--- Drug 1 vs Drug 3 ---\n")
## --- Drug 1 vs Drug 3 ---
# Create new survival object for this subset
surv_obj_1vs3 <- Surv(time = Traviata_1vs3$time, event = Traviata_1vs3$event)
# Perform log-rank test
drug1_vs_drug3 <- survdiff(surv_obj_1vs3 ~ group, data = Traviata_1vs3)
print(drug1_vs_drug3)
## Call:
## survdiff(formula = surv_obj_1vs3 ~ group, data = Traviata_1vs3)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## group=Drug 1 38 24 20.8 0.483 0.769
## group=Drug 3 45 34 37.2 0.271 0.769
##
## Chisq= 0.8 on 1 degrees of freedom, p= 0.4
p_value_1vs3 <- 1 - pchisq(drug1_vs_drug3$chisq, df = 1)
cat("P-value:", format.pval(p_value_1vs3, digits = 6), "\n")
## P-value: 0.380643
cat("Significant at α =", alpha_corrected, ":", p_value_1vs3 < alpha_corrected, "\n\n")
## Significant at α = 0.01666667 : FALSE
# 3. Drug 2 vs Drug 3
cat("--- Drug 2 vs Drug 3 ---\n")
## --- Drug 2 vs Drug 3 ---
# Create new survival object for this subset
surv_obj_2vs3 <- Surv(time = Traviata_2vs3$time, event = Traviata_2vs3$event)
# Perform log-rank test
drug2_vs_drug3 <- survdiff(surv_obj_2vs3 ~ group, data = Traviata_2vs3)
print(drug2_vs_drug3)
## Call:
## survdiff(formula = surv_obj_2vs3 ~ group, data = Traviata_2vs3)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## group=Drug 2 54 25 39.7 5.42 16.9
## group=Drug 3 45 34 19.3 11.11 16.9
##
## Chisq= 16.9 on 1 degrees of freedom, p= 4e-05
p_value_2vs3 <- 1 - pchisq(drug2_vs_drug3$chisq, df = 1)
cat("P-value:", format.pval(p_value_2vs3, digits = 6), "\n")
## P-value: 4.03402e-05
cat("Significant at α =", alpha_corrected, ":", p_value_2vs3 < alpha_corrected, "\n\n")
## Significant at α = 0.01666667 : TRUE
# Summary table of pairwise comparisons
cat("--- Summary of Pairwise Comparisons ---\n")
## --- Summary of Pairwise Comparisons ---
comparison <- c("Drug 1 vs Drug 2", "Drug 1 vs Drug 3", "Drug 2 vs Drug 3")
chisq <- c(drug1_vs_drug2$chisq, drug1_vs_drug3$chisq, drug2_vs_drug3$chisq)
pvalue <- c(p_value_1vs2, p_value_1vs3, p_value_2vs3)
significant <- pvalue < alpha_corrected
results_table <- data.frame(Comparison = comparison,
ChiSquare = round(chisq, 2),
P_Value = format.pval(pvalue, digits = 6),
Significant = significant)
print(results_table)
## Comparison ChiSquare P_Value Significant
## 1 Drug 1 vs Drug 2 26.08 3.27492e-07 TRUE
## 2 Drug 1 vs Drug 3 0.77 0.380643 FALSE
## 3 Drug 2 vs Drug 3 16.86 4.03402e-05 TRUE
# Conclusion
cat("\nConclusion:\n")
##
## Conclusion:
if(all(significant)) {
cat("All pairwise comparisons show significant differences in survival.\n")
} else if(!any(significant)) {
cat("None of the pairwise comparisons show significant differences in survival.\n")
} else {
sig_pairs <- comparison[significant]
cat("The following comparisons show significant differences in survival:\n")
for(i in sig_pairs) {
cat("- ", i, "\n")
}
}
## The following comparisons show significant differences in survival:
## - Drug 1 vs Drug 2
## - Drug 2 vs Drug 3